{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pylab as plt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "data = np.loadtxt('data100m.csv',delimiter=',')\n",
    "x = data[:,0][:,None]\n",
    "t = data[:,1]\n",
    "x = (x - 1896)/4\n",
    "sigma_sq = 0.05"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "testx = np.arange(0,40,0.1)[:,None]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x108373f10>]"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX0AAAEACAYAAABfxaZOAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAE05JREFUeJzt3XuMXGd9xvHnWSJQB1RqLrkoJg4INmkt4t1sa1GFOtOi\nEoOiJnEjSPgHshVyWwix1FYJbSU7UpEIfyBHINSQ2lZaYYUSC0hSUZwIJpaLIPZ4l1wax5FoTBKS\nJSoJF22Fop1f/5izy3izszvXnXPO+/1II8+8e87Me/bYzxy/5704IgQASMPYqCsAAFg/hD4AJITQ\nB4CEEPoAkBBCHwASQugDQELWDH3b+2zP2X6kpexztp+wPWv7kO3fbrPvdtsnbZ+yffMgKw4A6F4n\nV/oHJF2xrOywpM0RMSHpKUmfXr6T7TFJX8z23SzpetsX91ddAEA/1gz9iDgq6aVlZQ9GRCN7+X1J\nG1fYdaukpyLidES8IuluSVf1WV8AQB8G0aY/LelbK5SfL+mZltfPZmUAgBHpK/Rt/4OkVyLi4IDq\nAwAYorN63dH2xyR9UNKftNnkOUkXtLzemJW1ez8mAQKALkWEu9m+0yt9Z4/mC3u7pL+T9GcR8es2\n+xyT9E7bm2y/VtJ1ku5d7UMWJN04MaGFhQVFRGkeu3fvHnkdOD6Oj+Mr36MXnXTZPCjpe5LGbf/Y\n9g2SviDpDZIesH3C9peybc+zfb8kRcSCpE+q2dPncUl3R8QTq33WTVu2aOf+/RobY/gAAAzDms07\nEfGRFYoPtNn2eUlXtrz+T0kXdVqZ20+cIPABYIhylbBlDfxqtTrqKgwVx1dsHF9a3Gu70KDZjrzU\nBQCKwLZiSDdyAQAlQOgDQEIIfQBICKEPAAkh9AEgIYQ+ACSE0AeAhBD6AJAQQh8AEkLoA0BCCH0A\nSAihDwAJIfQBICGEPgAkhNAHgIQQ+gCQEEIfABJC6ANAQgh9AEgIoQ8ACSH0ASAhhD4AJITQB4CE\nEPoAkBBCHwASQugDQEIIfQBICKEPAAkh9AEgIYQ+ACSE0AeAhBD6AJAQQh8AEkLoA0BC1gx92/ts\nz9l+pKXsWtuP2V6wfekq+z5t+4e2Z2w/PKhKAwB608mV/gFJVywre1TSNZIeWmPfhqRqRExGxNYe\n6gcAGKCz1togIo7a3rSs7ElJsu01drdoQgKA3Bh2IIekB2wfs/3xIX8WAGANa17p9+myiHje9lvV\nDP8nIuJou4337Nmz9LxaraparQ65egBQHLVaTbVara/3cESsvVGzeee+iLhkWfl3Jf1NRJzo4D12\nS/plRHy+zc+jk7oAAJpsKyLWamY/Q6fNO84e7X62UmUqtt+QPX+9pPdLeqybygEABquTLpsHJX1P\n0rjtH9u+wfbVtp+R9B5J99v+Vrbtebbvz3Y9R9JR2zOSvq/m/xQOD+cwAACd6Kh5Zz3QvAMA3Rlm\n8w4AoAQIfQBICKEPAAkh9AEgIYQ+ACSE0AeAhBD6AJAQQh8AEkLoA0BCCH0ASAihDwAJIfQBICGE\nPgAkhNAHgIQQ+gCQEEIfABJC6ANAQgh9AEjIWaOuwHpqNBqamZmRJE1OTmpsjO88AGlJJvUen5nR\nrqkpnd62Tae3bdOuqSk9nn0BAEAqklgYvdFoaNfUlPbOzi59yzUk7ZqY0N56nSt+AIXEwuhtzMzM\nqHrq1BkHOybp8lOnlpp7ACAFSYQ+AKApidCfnJxUbXxcjZayhqSHxsc1OTk5qmoBwLpLok1fat7I\nvWN6WpefOiVJqr3rXfrLAwe0mdAHUFC9tOknE/oSXTYBlAuhDwAJofcOAGBVhD4AJITQB4CEEPoA\nkBBCHwASQugDQEIIfQBICKEPAAlZM/Rt77M9Z/uRlrJrbT9me8H2pavsu932SdunbN88qEoDAHrT\nyZX+AUlXLCt7VNI1kh5qt5PtMUlfzPbdLOl62xf3WE9JzWkU6vW66vW6Go3GwLYFgFSsGfoRcVTS\nS8vKnoyIpyStNvx3q6SnIuJ0RLwi6W5JV/Va0W5WvmKVLABY2TDb9M+X9EzL62ezsq41Gg3dMT2t\nvbOz2jE/rx3z89o7O6s7pqdfdRXfzbYAkJpC3MjtZuUrVskCgPbOGuJ7PyfpgpbXG7Oytvbs2bP0\nvFqtqlqtDqNeAFBItVpNtVqtr/foaGpl2xdKui8i3r2s/LuS/jYi6ivs8xpJT0p6n6TnJT0s6fqI\neKLNZ7SdWrmbhc1ZBB1AKoYyn77tg5Kqkt4saU7SbjVv7H5B0lskvSxpNiI+YPs8SXdGxJXZvtsl\n3a5mC8u+iPjsKp+z6nz63ax8xSpZAFJQ+kVUuln5ilWyAJRd6UMfAPAbrJwFAFgVoQ8ACSH0ASAh\nhD4AJITQB4CEEPoAkBBCHwASQugDQEIIfQBIyDBn2SyMvEzvMKz3ZkoKAIuS/9eflxW5hvXerCIG\n4AwRkYtHsyrra2FhIW6cmIgFKSJ7LEjNsoWFnrcdZj3y8L4A8iHLza6yNukr/bysyDWs92YVMQDL\nJR36AJCapEN/cnJStfFxtS6X3pD00Pi4JpctuNLNtsOsRx7eF0BxJT+ffl5W5BrWe7OKGFBeLKLS\nI7psAigiQh8AEsLKWQCAVRH6AJAQQh8AEkLoA0BCCH0ASAihDwAJIfQBICGEPgAkhNAHgISwclZB\nMbUCgF6QFAXEalgAesXcOwXTaDS0a2pKe2dnl76xG5J2TUxob73OFT+QEObeSQCrYQHoB6EPAAkh\n9AuG1bAA9IM2/QJiNSwAEouoJIUumwCGEvq290m6UtJcRFySlW2Q9FVJmyQ9LelDEfHzFfZ9WtLP\n1WyBeCUitq7yOYQ+AHRhWL13Dki6YlnZLZIejIiLJH1H0qfb7NuQVI2IydUCH8ijRqOher2uer2u\nRqOx9g5AAawZ+hFxVNJLy4qvknRX9vwuSVe32d2dfAaQNwyAQ1l11KZve5Ok+1qad34WEW9q+fkZ\nr1vKfyTpZUkLkr4cEXeu8hk07yAXGACHohjl4Kx2aX1ZRFwq6YOSPmH7vQP6PGBoGACHMut1wrU5\n2+dExJztcyX9dKWNIuL57M8XbX9d0lZJR9u96Z49e5aeV6tVVavVHqsHAOVTq9VUq9X6eo9Om3cu\nVLN5593Z69sk/SwibrN9s6QNEXHLsn0qksYi4le2Xy/psKRbI+Jwm8+geQe5QPMOimJYXTYPSqpK\nerOkOUm7JX1D0tckvU3SaTW7bL5s+zxJd0bElbbfLunrajb9nCXpKxHx2VU+h9BHbjAADkXA4Cxg\ngBgAh7wj9IE1EOQoE6ZWBlZB33uAK30kgpuzKCOu9IE26HsPNBH6AJAQQh9JYPEZoIk2fSSDvvco\nG7psAmugyybKhNAHgIT0Evq9TriGkurmSjgvV815qQdQBPzrwJJuBi/lZaBTXuoBFAXNO5DU3eCl\nvAx0yks9gFFhcBZ61s3gpbwMdMpLPYAiIfQBICGEPiR1N3gpLwOd8lKPFDQaDdXrddXrdTUajbV3\nQG7Rpo8l3QxeystAp7zUo8wWf8fVxd/x+Lh27t/P7zgH6KePvtFlE624WZ5vhD6AgarX6zq9bZt2\nzM+fUX6oUtGFR45oampqRDWDxOAsYGT43waKgr+ZQJ/KPECMm+XlQ/MO0IcU2ry5WZ5ftOkD6yyV\nNm+ar/KJNn0AQzE2NlaaL7DU8XUN9IE2bxQNzTtAn/LU5k0zTFpo0wdGJA9hy8jZ9BD6QKJS6EWE\nV2NqZSBRTDONThH6AJAQQh8oAXoRoVO06QMlkadeRFgf3MgFElfEqbHRO0IfQEfo3lkOhD6ANdG9\nszzosglgTXTvTBuhDwAJWTP0be+zPWf7kZayDbYP237S9rdtv7HNvtttn7R9yvbNg6w4UFSNRkP1\nel31el2NRmPtHQaM7p1p6+RK/4CkK5aV3SLpwYi4SNJ3JH16+U62xyR9Mdt3s6TrbV/cX3WBYsvD\nKltjY2PauX+/dk1M6FClokOVim7askU79+8fSHv+qL/UsLqObuTa3iTpvoi4JHt9UtLlETFn+1xJ\ntYi4eNk+75G0OyI+kL2+RVJExG1tPoMbuSi1vN1AHUaXTXoFra/1vJF7dkTMSVJEvCDp7BW2OV/S\nMy2vn83KgCTl7Qbq4sIoU1NTA7vCv2N6WntnZ7Vjfl475ue1d3ZWd0xPc8WfI4NaOWsgl+h79uxZ\nel6tVlWtVgfxtgDWwVpfav2uvMVgMqlWq6lWq/X1Hr2G/pztc1qad366wjbPSbqg5fXGrKyt1tAH\nymZyclJ3jY/r6mXNOw+Nj+uaNs0fBF3T8majuxJtNlp+MXzrrbd2/R6d/g1y9lh0r6SPZc8/Kumb\nK+xzTNI7bW+y/VpJ12X7AUnq9gZqHm76dmNYvYJoNhqwiFj1IemgpJ9I+rWkH0u6QdIGSQ9KelLS\nYUm/k217nqT7W/bdnm3zlKRb1vicAFKwsLAQx48fj+PHj8fCwkLbbW6cmIgFKSJ7LEjNsjb75MFj\nJ07EjRMTcU+lEvdUKvHJLVvisRMn2m7fye/i+PHjcahSWfo9LD7uqVTi+PHjwzqUQshyc80cb30w\nDQOQQ/V6Xae3bdOO+fkzyg9VKrrwyJG+28eHqdMmqU57+hT5dzFsTMMAYOQ66RXUTZMNg8kGi9AH\ncqjsQddN99VhDyZLzaC6bAIYoKWgW74oSqJBt3lyUnvr9aUvhNsT7snUL9r0gRwra5fNFEYnrwfm\n0wdQGHlZ3rHIU0cQ+gAKZZhX2J28dy//48jT/wrovQOgUAY9/8+iTge2dTsfUtEGzK2EG7kASqW1\nO+himF89O6td09N93S8Y1vuut2LUEgA61M3VezddY/M2S2qvuNIHkKwUu8ZyIxdAqQzr5mzeuplK\n9N4BAEnD6w6al26miwh9AMgMq2tl0btsEvoAUFD00wcArIrQB4CEEPoAkBBCHwASQugDQEIYkQsA\nQ5Kn7p2LRl8DACihvM7IST99ABiw9ZqygX76AJADeZ6Rk9AHgIQQ+gAwYN3M07/080ZD9Xpd9Xpd\njUZjxW0GgdAHgAFbmqd/YkKHKhUdqlR005Yt2tlmnv71vOnLjVwAGJJhz9PPLJsAUDD1el2nt23T\njvn5M8oPVSq68MgRTU1Ntd2X3jsAgFUR+gAwQr3c9O0HzTsAMGK9LsNImz4AFFQv8/QQ+gCQEG7k\nAgBWRegDQEL6Cn3bN9l+NHt8aoWfX277Zdsnssc/9vN5AID+9Bz6tjdL+gtJvy9pQtKVtt+xwqZH\nIuLS7PFPvX5ekdVqtVFXYag4vmLj+NLSz5X+70r6QUT8OiIWJB2RtGOF7bq6yVBGZf9Lx/EVG8eX\nln5C/zFJf2R7g+2KpA9KetsK2/2h7Vnb/2H79/r4PABAn3peIzciTtq+TdIDkn4laUbSwrLN6pIu\niIh52x+Q9A1J471+JgCgPwPrp2/7M5KeiYh/XmWb/5E0FRE/W+FndNIHgC5120+/5yt9SbL91oh4\n0fYFkq6R9J5lPz8nIuay51vV/JJ5VeBL3VccANC9vkJf0iHbb5L0iqS/johf2N4pKSLiy5Kutf1X\n2c//T9KH+/w8AEAfcjMNAwBg+EY+Itf2dtsnbZ+yffOo6zNotp+2/UPbM7YfHnV9+mV7n+0524+0\nlG2wfdj2k7a/bfuNo6xjP9oc327bz7YMMtw+yjr2yvZG29+x/XjrgMqynL8Vju/GrLws5+91tn+Q\nZcmjtndn5V2dv5Fe6dsek3RK0vsk/UTSMUnXRcTJkVVqwGz/SM2b1y+Nui6DYPu9avbW+teIuCQr\nu03S/0bE57Iv7g0Rccso69mrNse3W9IvI+LzI61cn2yfK+nciJi1/QY1e9ddJekGleD8rXJ8H1YJ\nzp8k2a5kvSFfI+m/JH1K0p+ri/M36iv9rZKeiojTEfGKpLvVPEllYo3+9zwwEXFU0vIvsKsk3ZU9\nv0vS1etaqQFqc3xSCQYZRsQLETGbPf+VpCckbVRJzl+b4zs/+3Hhz58kRcTimoqvU/OebKjL8zfq\nMDpf0jMtr5/Vb05SWYSkB2wfs/3xUVdmSM5e7KUVES9IOnvE9RmGT2aDDP+lqM0frWxfqOb0Kd+X\ndE7Zzl/L8f0gKyrF+bM9ZntG0guSHoiIY+ry/I069FNwWURcquaI5U9kzQdlV7beAV+S9I6ImFDz\nH1uhmwmypo97JN2UXREvP1+FPn8rHF9pzl9ENCJiUs3/oW3N5kDr6vyNOvSfk3RBy+uNWVlpRMTz\n2Z8vSvq6mk1aZTNn+xxpqV31pyOuz0BFxIstK/zcKekPRlmfftg+S81A/LeI+GZWXJrzt9Lxlen8\nLYqIX0iqSdquLs/fqEP/mKR32t5k+7WSrpN074jrNDC2K9lVh2y/XtL71ZyzqOisM9tI75X0sez5\nRyV9c/kOBXPG8WX/kBbtULHP4X5J/x0Rt7eUlen8ver4ynL+bL9lsWnK9m9J+lM171t0df5G3k8/\n6z51u5pfQPsi4rMjrdAA2X67mlf3oeZNl68U/fhsH5RUlfRmSXOSdqs5p9LX1Jxw77SkD0XEy6Oq\nYz/aHN8fq9k+3JD0tKSdi22oRWL7MjVnw31Uzb+TIenvJT0s6d9V8PO3yvF9ROU4f+9W80btWPb4\nakR8Jhsg2/H5G3noAwDWz6ibdwAA64jQB4CEEPoAkBBCHwASQugDQEIIfQBICKEPAAkh9AEgIf8P\nstuG+ZPfl7oAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10599e090>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure()\n",
    "plt.plot(x,t,'ro')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define the prior mean and covariance for $w$ ($k$ is the polynomial order)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 0.]\n",
      " [ 0.]\n",
      " [ 0.]\n",
      " [ 0.]\n",
      " [ 0.]]\n",
      "[[ 100.    0.    0.    0.    0.]\n",
      " [   0.  100.    0.    0.    0.]\n",
      " [   0.    0.  100.    0.    0.]\n",
      " [   0.    0.    0.  100.    0.]\n",
      " [   0.    0.    0.    0.  100.]]\n"
     ]
    }
   ],
   "source": [
    "k = 4 # linear\n",
    "prior_mean = np.zeros((k+1,1))\n",
    "prior_covariance = np.zeros((k+1,k+1))\n",
    "for i in range(k+1):\n",
    "    prior_covariance[i,i] = 100\n",
    "print prior_mean\n",
    "print prior_covariance"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create $\\mathbf{X}$ and $\\mathbf{X}_{test}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X = np.ones_like(x)\n",
    "Xtest = np.ones_like(testx)\n",
    "for i in range(k):\n",
    "    X = np.hstack((X,x**(i+1)))\n",
    "    Xtest = np.hstack((Xtest,testx**(i+1)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sample some $\\mathbf{w}$ from the prior and plot the functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAEGCAYAAABxfL6kAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXl8VOXVx7/PLEkIWVnCFlYJArIjIKAQUBarosW9Vi0o\nal2rfdvqqxax9u1r9VVrq23dWrUKiiu4IQhxQRFUdtlkJ2whEJbsM/O8f5wZEjCBZDJbkvP9fO7n\nzsy9uffMZfjdc885z3mMtRZFURSl4eCItgGKoihKaFFhVxRFaWCosCuKojQwVNgVRVEaGCrsiqIo\nDQwVdkVRlAZGRIXdGPO8MWaPMWZFDfZ9zBiz1BjznTFmnTFmfyRsVBRFqe+YSNaxG2POBI4AL1lr\n+9Ti724F+llrrw+bcYqiKA2EiHrs1tovgAOVPzPGdDHGfGiMWWKM+dQY062KP70SmB4RIxVFUeo5\nrmgbADwD3Git3WiMGQz8HTg7sNEY0wHoBMyPjnmKoij1i6gKuzGmKTAMmGmMMf6P3cftdgXwhtXe\nB4qiKDUi2h67AzhgrR1wgn2uAG6OkD2Koij1njrH2I0xmcaY+caY1caYlcaY20/2J/4Fa+1hYLMx\n5pJKx+tT6XV3IM1au6iudiqKojQWQpE89QB3WWtPA4YCt/gF+UcYY14FvgS6GWO2GWMmAVcB1xlj\nlhljVgETKv3J5cCMENioKIrSaAh5uaMx5h3gr9baT0J6YEVRFKVGhLTc0RjTCegHfB3K4yqKoig1\nJ2TCboxJAt4A7rDWHgnVcRVFUZTaEZKqGGOMCxH1l62171azj5YrKoqiBIG11px8rwpC5bG/AHxv\nrf3LiXay1sb8MnXq1KjboHaqjWqn2hlYgiEU5Y7DkcqW0ZWado2v63EVRVEaM2vXws6dwf1tnUMx\n1tqFgLOux1EURVHAWnjmGbjvPvj3v4M7RrRHnsYc2dnZ0TahRqidoaM+2AhqZ6iJRTt374YpU2Db\nNjj3XCguDu44EWvba4yxkTqXoihKfeONN+DWW+EnP4GFC2HgQHjqKWjWzGBrmTxVYVcURYkiBw7A\nbbfB11/DiBHw3nvw5JNw+eWy3ZjaC7tOjacoihIl5s6FPn0krp6cDLt2wdKlFaIeLCrsiqIoEebI\nEbjlFpg8Gc4/H+bMgRtvhPffh7Zt6358FXZFUZQIMmcO9OolpYyZmbBiBSxaJMJuahVwqR4VdkVR\nlAiQnw/XXCMCPno0fP45/Oxn8Nln0LVraM+lwq4oihJGrIXXXhMv3eOB1FTIzYVvvpGkqTMMo4C0\nKkZRFCVM5ObCzTfD+vVSvjh3LjzyCFx9dc3DLloVoyiKEgN4PPD449C3r1S7eDxQXi7x9GuuCV0s\nvTrUY1cURQkhX3whXnpqqoj6mjVSl37BBcEdTz12RVGUKLF3L0yaJDXo/fqJoA8cCKtXBy/qwaLC\nriiKUge8XvjHPyQ5WlQEaWnS8+XLL+EPf4DExMjbpKEYRVGUIMnJgTvvBLcbWreG776T2Poll4Qu\njq6hGEVRlAiwcSNMnCiJ0G7d5P1pp0n45dJLw58cPRkq7IqiKDXk4EH4zW9g8GCIi5PPrJWa9D/9\nSZKlsYAKu6IoyknweCSOfuqpMrNR586wYQO88gq8/rq8jyU0xq4oilIN1sKbb8K990pSND0dVq6E\nP/5RwjCOCLjGwcTYdQYlRVGUKvjkE7j7bigpgawsadR1550yIUZSUrStOzEailEURanEt9/CmDEy\nRV27dtKFsXt3WLdOPPdYF3VQj11RFAWAVatg2jTpujhsmPRMb94cli2D9u2jbV3tUGFXFKVRs3Il\nPPigtM8dNEji6g4HfPop9OgRbeuCQ0MxiqI0SlaskIFE55wDxcXyWZMm0oHxjTfqr6hDiITdGPO8\nMWaPMWZFKI6nKIoSLpYtE0EfM0YSoyDD/ufNg5kzZQ7S+k6oPPZ/AeNCdCxFUZSQYi3Mnw/jxsG5\n58Lhw/JZ06ZS/fL669C7d7StDB0hibFba78wxnQMxbEURVFChdcLb70FDz8MBQXQqROUlkoy9Mkn\nZcBRQ0STp4qiNDiKi+HFF+HRR8UrT0+HzZvh4ovhpZegbdtoWxheIirsDzzwwNHX2dnZZGdnR/L0\niqI0cHbulKH///yneOfNm8P27TJp9LvvyuQXsU5OTg45OTl1OkbIWgr4QzGzrbVVph60pYCiKOHA\nWvjqK/jrX+Gjj6Qv+vbt0KwZ3HGHTHyRkBBtK2tPWZksycnRbSlg/IuiKErYKS2F116TWHleHnT0\nZ/natJFOi8OHR799bm3xeqWe/tVXJTfw+OPBHSckwm6MeRXIBpobY7YBU621/wrFsRVFUSqzaRM8\n95ws7drJJBdFRSLkr7xS/0aJ+nzSh2bmTKnOyciQ0NHSpdChA1x7be2PGaqqmJ+F4jiKoihVUVoK\n77wDzz4rsxR16yajQxMSpKfLlVfK4KL6gs8nU+fNnCndI1NTpbZ+7lzo2bPux9eqGEVRYpa1a0XM\nX3pJQixOp4RXBg2CZ56pX4OJvF744gsZ1frmm9CihYj5xx+HRswro8KuKEpMUVAg4vfiiyLsXbuK\nd56aKt75JZdEZ4LoYCgthQULYNYsiZm3aSP2L1gQ3hp6FXZFUaJOWZlUtLz8sqy7dIHCQpl+7qyz\n4PnnQ+/Vhou8PPjgAxHzefOkSmfCBOkamZUVGRtU2BVFiQrWwuLFIuYzZkhoIi5Owi19+8oMRaNG\nyftYxlqZxHr2bFlWrpTGYhMmSE19y5aRt0mFXVGUiGGtNOF6/XVZiouhVSuZU7RDB7j6avjpT2N/\nMosjRyAnB+bMgQ8/lCeOCy6A++6D7Ozo182rsCuKElaOF/OiIvHO9+0Tz/yyy2Sof5s20ba0enw+\nWL5chPzjj2HJEkngjhsn+YC+fWOrZl4ns1YUJeT4fCJ+77wjg4iKiiQksX27dFG89FIR83btom1p\n9ezeLTHyOXOkDDE1FcaOFTHPzo7cU0Uwk1mrsCuKEhKKikQIZ82Svixut4jfzp1w+ulw0UVSEZKZ\nGW1Lq2bvXgmv5ORI1cqePSLg48aJoHfuHB27VNgVRYkou3bBe++JmM+fD61bi7eenw/jx0sC8dxz\npRlXrLFvn0x/t2CBiHlurlTgZGdL0rZPn9hI3KqwK4oSVsrKZMTkxx9L0nDDBmmBe+gQuFzilU+Y\nACNHQnx8tK2twFppRbBwYcWyfbu0IRg1Spb+/WNDyI9HhV1RlJBiLaxbJ0I+Z454tunpIoB5eTB0\nqIQqxo0TDzdWEoilpdJrJSDiX34pN57hw2HYMFn37SvholhHhV1RlDqzbZuEKHJyxCsvLhYx37dP\nKld+8hOJOY8cGRsjQK2FrVslWbtkibTwXbpUBgMNH14h5h06xM6NpzaosCuKUiushS1bKoR87lwJ\nq6SlyToxUQbbZGeLmMdC58TduytEfMkS+OYb8bwHDZLljDNgyBBITo62pcGRnw+rV8syZAgMHBjd\nfuyKosQ4Ho8IxldfyRD3Tz6RapbUVDh4UMRw4kQR8pEjpRIkWl6utbBjh9SPL18O334rQl5UJFU2\ngwbBTTfJ61gum6wKa+UJaO1aGbW6ejWsWAGrVkkrhRYtJEfxy18Gd3z12BWlAbNnj/T6XrRIvPKl\nS8ULd7ul2VaHDlIJMmKECHmnTtGxs7QUvv++QsSXLROhc7slFt63LwwYIGLepUv9CamUl0vSdu1a\nyVWsWiXfa+NG6faYnCzfpbhYlk6doEcPaUvcvTuMHg1dumgoRlEaLQcPiiB+950kC7/4QsQ7LU1E\nw+eTMMVZZ0nSc9CgyM8B6vVK6Of778VTXbFChPyHH+CUUypEvF8/WbdqFVn7gsHnk1r9jRvleyxf\nLt9rwwapjQ/cSEtLoaREqoi6d5emZt26SS4gK0vCXA7Hj4+vMXZFaSTk5YmAL10qYZUlS+TRPiVF\nwi3FxeL5DR0qyxlniHhEytMtKxORW7NGRDwg5OvXywxBPXqIsPXqJQLes2f0+6uciLIyuSFt3Chh\nk+XLxQPftk1i4k6nLGVlsm7TRp4sevYUEQ+Id4cOtS+pVGFXlAZGcbE8xgfEZPHiijhsUpIIidcL\np50mlR+nny4hi27dpLwvnPh8EgP/4QfxTn/4QZZ162DzZvFAe/asEPGePaUHeSw2+CorkwFKmzbJ\ndV6zRkR8+3YJZxUWyvX0+WT/5s1lBG1WlpR5nnqq5CM6dZIKolCiwq4o9ZSSEhHFQJx58WIR87y8\nipLCoiIJTfTuDYMHy7p/fxGUqh7hQ4HHI4JXWbgDrzdtgmbNZCKMrKyKdVaW3FhiyQM/dEiEeuVK\nWTZskBLJPXskXFVSItfQ5xO709PlWnfsKDemXr3k+3XuLInNSMb4VdgVJYbxekVM1q8XL3zpUvG+\nt2yR+HhCQkUirXlz8cIHDZJQxWmniVcYarEsKxOvdOtWsaPyeutWEfWMjB+Ld9euEhNv2jS09gRj\n/8aNch3XrpWbzbZtItj5+SLoJSUVnnZcnOQVWraUJ4pAuKR3bxHtNm1ib9CSCruiRJmCAhHGzZtF\nwFetkvW2beJ9x8eLZ1hSIp54hw4iLAMGiHB36yZiEwoBLyuTXi65uT9etm0TO/fulWRex44SRjh+\nnZkZ+dYAXq8I9Pr1ItqBG8zu3ZJH2L+/QrC9Xvmb+Hi5yaSliWi3ayf2d+smHndWltygYrFlwMlQ\nYVeUMOLxiLjs3Cmx5TVrjo3F5uWJZ+h2i+B4vRKqaNdOPNzevUVoApUQwQ6gKSmRc+3ZI8JclXjv\n3AkHDoiYtWt37FJZyNu1C28s3loR4cCNZMcOuVa5uWJ3Xp7YefiwhJrKyuQaGiN2NWkiMfm0NHmK\nadVKboZZWRXJ1/T0+lP+GAwq7IoSBD6fiEtAIDdtkhjspk0iRHv2yPaiIhEbY0Tk3e6KWGzXrhUV\nEIEkWqtWNYt9ezxy/Pz8YwV7796K15XXxcUi2BkZco7WrSsEu7KAh8pDtVZmDDp4UM6/Y4csAS86\nL09sLygQgS4slJtPILELFUIdHy9PKqmpctPLyKi40XTuLIIdqwnWmmCtlDUWFcl1KCqSa3f48LHr\nqj6rav3ggzBpUpSE3RgzHngCcADPW2sfrmIf6/V6cRz3S/f5fCxduhSA/v37/2h7VfsAJ3wfzDlO\nRjB21taOUJwjFETCznCew+sVkdm3z8cXX6xnzx43xnRmxw7HUTHKz4eDBy2HDllKSw1Op3hGAW8x\nKUmEp3VrS3LyfjIzyxg5shUdOzpo21bEqGnTChushe7d+3P4sIODByVcsH+/2LByZS4HD7qIi2vN\ngQOG/PyK7fv3w5EjlqQkL6mpHtq3jycjw9CqVYVwt2jh4+DBDTRrVk52dk+aNXP8yEOt6lp4vSIu\nhYVw8KCPhQvXk5/vJjGxM/v2OcjLk/MfOCDXYu/eEoqLHXi9cZSUGEpLZYBNID4NFmPA6bQkJBgS\nEw1JSVJi2ayZJBVbteLo9QmIdatWsdHp0Vq52ZSUyM2xpORYAT7Rujb7uFxy82raVNZJSfJ0lpR0\n7OuarFu2hKZNoyDsxhgHsB44G9gJLAGusNauPW4/e1u/ftz4wguc5hfj1UuX8s/Jk8levx6AnG7d\njtle1T7vZGaSAIzfsaPK98cfoybnOBnB2FlbO0JxjlAQCTtPtn35kqU8Pfl2Bv6whyKbzOeZfRk6\n5X7cCZ2PerKbfihg7ZKNuIriKLYpHHamYNwplHuceDzgcFisz4fF4MDiMGWkpLto0cJNy5bQJL6A\nPcs/pPfB70k3u9ncLokJ/3UL3Xt1xeeT8MH3K7bw0T9nkLm3mCKbyg+p7ek4+GxwNOPQIdizu5jd\nW/fjK0ugmBSM8ZDazEXz5m6aNQO36yB7V39O18NbSTL72Z3h4rxbrqTPgE6kp4sgbt+4khm/vZNB\nm3IptQks7NCTMb96gLSWWRQUwLrVuXw64wNa7iuniBS2N21Ly+6nYxwpR4Xk8OEyDuUXYrxOPLjx\n4ALjwlrjvwFYrAWDxYEXYzw0TXaTlOyiaVNwOQs5uH01mUU7STIHKGhuGXXVWPoOyKRtW0ko5u1c\nyWt3XcfZP6wO6rfn88lNorRUxLWs7NjXld+XllYIb2URruvrkpKK8E5CgiyVBbiq9Ym2VbdvqENb\nUQnFGGPOAKZaa8/1v78bsMd77cYY6wV+1a8fT3z7LQC/GjiQJ5YtI+CL+Sptdzgc+Hy+Y/bxAXcA\nf4Eq3x9/jJqc42Qcb0NN7KytHcFci2C+S12/K8AdAwbyyPLVlNKEYppwhCbc320I1z/9KocPOzhw\nwMe/7/8zY3IPUUgSR/zLF6ntOPXMsRw5AmsWr6ZlsaWYRIppQgkJHCIB42qCxxOwJvC7NIDF4CMu\n3onbbYiLs5QcyqOV5xBNKPFbUszWJBd9Rg7F7YYVCxbQ8+ARPLgppgmFJLIxvgVJGZ0pLIRDB0rx\n2jji/d8knhIK3dC6UxtcLkN8vGXnulV0LT6IEy9OvBh8bEhK4tQhgwFYt2gxXQsL8eGkHBflxLE5\nIY3mHbIoK4O8HXtJ9BjKcVOOG49/PxwufD7j/44Wc/QbG/839uFyO3C5wFNaSLKvhDjKSaCEBErY\nlwDdBnQnMdGQmOhjzWdzGFGwjWQKacoRmnKEOe3SuObBu4mLg//c9zuu3/oD4MSDi1LcvNy+C5fe\nOxWPB17/06NckLsHLy48uCjDzQet2pF99fWUlxtKSiyfz3yX/vsPU04cpcRTShyrkprTuf9gSkvN\nScXa45GKlLg48d4Dr6t6HxDdhIRjRbiur+PjNXlam5NeDIyz1t7gf/9zYLC19vbj9rPDmU8eTpr0\n7SuxqBXLaIEXizn6k96Hk/g+fUhKSuXIkQJKVqykBV4McBhDCdAc+U9QBJRgaUbgv78DB5CHIb5X\nHwDKVq2kGV7/HkI+DuJO60VS0o/HUxsjj2wg68LCQ5SuXkUzfEfPYjHsxxDfsxeJiSkUFR2i9PvV\nNMceI0f5OInr0QOA0jVr/XYcPRP7cRLXvTvWQsm6DaRj8WH8MuagACfuTl2Ij0+kqKiI8u25JPq/\np9fvfxXixNGiJS5XAtaKzT5fzdY/XnzYo8/dppLkBL7RSX8NBMTq+L+yx70zP/qs8icWc/QYlfZ1\nOAM/cqzPe/Q3Ezi+xfi3V5yxtnYbI8cACz4PTr8tgcWHxZEQL7+T4iLcR//FLA4s5VicqSkYA76C\nAyT4twe+RwkWZ/MWuN3xlJeX4cnfT5z/31sWJ6U4sfFN8XgMXq/xy7Esbsrx4SGuRSpNmsTh85Vi\nd24lzcrtI7BPgcNH0oDeMhpyyVe095Uc3ebCQ64T0i4Yh9MJxe/O5FRPydFtLjyscxsybppMly6Z\n7N69jUOP/y9Dyw4TR9nRZXG8IfOv/0fv3qeeVKzd7oad4AwnMS/s8PtKn2T7l2qPTGWBCA81uVYn\nsyFgZ03tPdF+lWXw+Nc1saOuVHWO4wXSVpLSirVxODDGgbW+SoLoqySIPlyJTSREcuQwCXhx4MWJ\nDwdeSoC41hk4HBbvrlxSrMcvc+ItHzYQ360rKSlNKSk5Qvnq5WT4KoTIiZc8hyVpxDCcTihdMJfO\nvhLclBHnF6zNLkP6FZficFiOvPoyvT2FxFN61J9e4XbR8q476Ny5Azt2bGH/nx9mcFnRURsc+FgS\n56L1Hx/E4bDk//dvObO0+Oh2J14+j4+j7VNP4nTC3pumcE5p4THfY15CHO1f/Q99+/bi++9XsuvS\niZxXUnjMMWY3SeCUue/jclm2jxrJJcWFx/yLvJmYSKfPPmPgwIF8++23bB0xgolFRVXuA9Rpe03O\nMXDgwBP/rJRakZOTQ05OztH306ZNq7Wwi/dThwU4A/io0vu7gd9VsZ/1gr2tXz/r9Xqt1+uV15Uc\nxsrbrbU/2scL9lb/uqr3xx+jJuc4GTW189a+/WwZDluK2xaRYAtIsjf0OtPu2eO1O3d67aTTxtjN\ntLcb6WTXkmWX0cte3u1yu3Ch1+bkeO1FXabY9xhn3+F8O4NL7fNcY0e0/aP905+8dupUa//rv7y2\nd/NX7HX8w/6cF+0lvGbP5x3bPmmhHTHCZ4cOtXbAAGt79bI2K8vaDh2sbd3a2vR0a5OSrI2Pt9bp\nlK/gdFobF2dtQoK1TZpY27Sp7JOc7LNuxyGbwgGbTIFN4qBN5LB1miIbF+ezbrfPQrkFjwWfNXit\nA481lFuXy2ddLmtdLp91UGZd/sVJuXXgseC14As8Uxy3eK3BY91un01I8Fm344hNJ99msMu2Z6vN\nYq1t0WSNHTbMZ8eOtfbii322R7O37K/4PzuN++yT3GJf5gp7YZcb7eLFXrt2rddOPu0cW0CS9WKC\n+l2cbJ+6bg/FOWLFTiW8iEzXTpdD4bE7gXVI8nQXsBi40lq75rj97K19+3LTv/71oyTayEASLSvr\nmO1V7fNuu3YkGMM4f7L0+PfHH6Mm5zgZwdhZWztCcY6a4PMdm8U/cqRiffAgrF6+jbkvvEmrveUU\nksqmpA60Pm0YHl8qBw7Avn3lHMj34rCQYvZBfCFderWlXWYyKSmSzS8+soeVc2bROX8vDjxsaZFO\n3/MuIz6xNYcOwfatB1j3zQbiC10U2WQOudLAnU5JqUse3d0ePMUHSPQU4zAeStxuUjIycLvj8fkk\nXltSUk5hwREcXoMXN+XGDcaNz2f8JYYWn88HOHDiwZhyklLiSE5xkZwMLucRCrYspUNhLsmOfPa3\n8HH+TRdw+pBOpKVJ3fS+nSuYcdckRm1YG7V/04b021OCI2p17P5yx0AO83lr7f9WsY+WOzagcsev\nvlrGgQNuMjJOO1o6F6hYkVpsy9atRezf72b/fjfp6YbMTKmvzsyEtm19eL1badmynJEju9Kxo4PE\nRKnd3b8f8vN9LF78A4cOuUhN7URBgeNoaV5FqaBl504PBQUuUlMhI8OQkcHRipPERB+HDuUBlvT0\nDMrKHJSUyA0sUOKXn19GcbGD8nIXxcVS7hioVff5wOu1uN2WhAQfzZo5adbM0LKllPDJd/FRXLyJ\ntDQPZ5zRjTZtHMcMmInE7yYUx4iV357yY3SAkhKTeL0i+IFBLcevc3NlNGLTpjKcvqolM7P6igaf\nT4S+ugE9gfWePTIiMyWl4gZT+WbTrp2MbmzSRKo59u2Tv922TWzdtUveSw28lEQWF8uNIFDiVl4u\n3zclRY7VurUcOzBMv3XrY5cmTSL2z6DUU1TYlXqLtSK8mzZJn5VNm45d8vKkadOpp0rvj8pLWlrN\nz+PziWBXvrlUXgI3GadTBtdUXrp0qRhVGmh+5fXK8XbvFuHftUtuBJs2yXF27ZLtBQWyf6DkzuuV\nm0J8fMUw+UD3wPbtK5bMTBX/xo4Ku9JgKSmRXiOBOSIDy9q1IrIBke/VS2bf6dMn+GHp1kqoZsuW\nihtN5WXrVvHIO3eWDofduskNJ9DEq6qOh9ZKmCnQK2XbNjnOxo0VvWb27auo4zZGvP/CQhn00q6d\nnOv4J5nOnaPfYVEJLyrsSqPDWvGyA0K/cqVMD7d6tXi7/fodu7RpU/d6ap9PPPHNm0WU16+XySXW\nrZM+5S1aVAh9YOneXbzyE53b55NQT2Xh37RJvtcPP8g5A6McAx0iAxNQB24wlUW/a9fQfF8luqiw\nK4ofj0eEdtmyYxeHA4YMkanihgyRfucpKaE7r9crorxuXcUExoHXhYXyRNG7tzxR9O4tS01DSV6v\nePwBLz8wx+a6dXIDcDolcRwXJ97+gQPyN927y9NM5ZtMVlZsTYShVI8Ku6KcAGtFGL/+GhYtkmXZ\nMglnBMT+zDNF/MLh5ebnV8zgs2KFrFetkgZaAbEfMECmt+vUqXY2WCtx/kCoKrD+/ntJLGdkiKcf\nCAnl5VVMqnzqqTKRR+/ecuOpr50VGyoq7IpSS8rLRWS//lomhf78cwlxjBwJ2dmy7tEjfOEMn09C\nOgGx//Zb+OYbseH0049dMjODs+Pw4Yqnhu+/l5vJihUS08/MrHhiOHBAnjbatKl4mgjccLKywj+H\nqlI1KuyKEgK2bIFPP5UlJ0dCKCNGwKhRMG6cxLPDza5dFSL/zTewZIl8fvrp8nQxbJisg52sA6RS\nJyDyy5fLetUqCU21aSPx/LIyeRLYu1c8+379YOBAWfr2rZiPVQkfKuyKEga2bRORnz8f5swRMRs/\nXpbs7MiELgJJ4iVLJIS0cKHMmXrqqTB8eMXSvn3dzhN4gli+HL77ruLGEhcnydi0NHnKyc2V+P4p\np1QI/cCBIvwq9qFFhV1Rwoy1Ejb56CNZliyBwYPhJz+Biy6KjDcfoLRUvPqFC+HLL2UdHw9nnSVP\nF6NHS3VMXcNI1kqFTkDkA0tqqlTipKeLZ791q4R8unSRazJ0qDxZ9OhRs5mklKpRYVeUCHP4MCxY\nAO+9B7NmyYw3P/2piHz//pEtNbRWvOjPPhOb5s+XdrmjR8syapTE1EOBzydVOV9/XXFj2bhRvnPX\nrlJxs2+fPFXs2ydho6FDZRkypHaDyho7KuyKEkW8XhG6d96Bt98Wj/qii+DSSyVMEmmv1VrxoAMi\nv2CBVOCMHi1hpLPPrluM/ngOHZLv/+WXsixaJG0TBgyQdVmZjC/49ltprzB8eEWCum3b0NnR0FBh\nV5QYwVopN3zrLXj9dUlUXn45XHll5D35AD6fhJHmzZMw0qJFFWGkc88NffWP1ytC/vnnclPJyRFP\nfeRICdd4vVJu+umn0lenciVSqJ4sGgIq7IoSo6xeDdOny+JyicBfeaUkP6PFkSPiyX/4IXzwgYh6\nQORHjw59qwKfT67DggWyfPqphK6ysyV84/VK7P7TTyV+P3KkhI/GjJF+Oo0VFXZFiXGslYTr9Okw\nY4YMjpo8Wbz5UIZFgrHr++9F5N9/X8IlZ58t+YLzz5cQTqjx+aTEcsEC+OQTyQ307Aljx8oNLz+/\nohqpY0cpNR07VkI4jWnUrAq7otQjPB4JiTz/vIQpLroIrrtOhCva/V3275eE8Ntvi7AOGlSRFG7X\nLjznLC3+/VAKAAAYL0lEQVSVROxHH0lZ6fbtcM45srRqJYnYjz+WWvszzxSRHzs2vAPIYgEVdkWp\np+zZAy+/LCLv9cL114vIN28ebctkgNacOSLy778vJY5XXAGXXRbepOfOnSLkc+bA3LkSd7/wQgnP\n7Nsnn82ZI9frvPPgggskhNTQ2hyrsCtKPcdaSWr+85/w7rswcSLcdpsM/IkFysvFg58+Xezr319y\nBRdfHJ5wTQCvV1o+zJol5y0sFCGfMEGeID7+GGbPlmTsqFGy7bzzpBqnvqPCrigNiLw8ePZZ+Pvf\npSnYbbdJOMTtjrZlQkmJJF2nTxdhHTFCRP7CC8PfI37dOhH5WbOk0ueccyRMNHy4hHNmzxabunUT\nkf/pT6XRWX1EhV1RGiAej9TG//WvMgDp1lvhl7+MrUE+hw6JJz19unjWEydKUnjYsPDHv/PyJET0\n1luSbB05UsJE48dLa4RZs2RbcjJccoksvXvXn7i8CruiNHCWL4dHHxVP+brr4Fe/ir3BPbt2Sb7g\nX/+SypdJk+CaayJj58GDIuSvvy4iP2qUiPx554mX/8YbMHOm9L4JiHy0xhXUFBV2RWkkbN0Kjz8O\nL70k3vFvfhPdmviqsFZGor7wggjq0KEi8hMmiLCGm4KCCpH/7DNJrF51lZRvrl5dIfLWyujgq66S\nFsWxhgq7ojQy8vPhqadkyc6GqVOlFjzWKCqCN9+Uqp/166Xq54YbIjfCtKBAqnr+8x8pm7z4YnmK\nGDZMYvQzZsCrr0p46+c/l1xBXTtlhopghF17rilKPaZ5c/j976UB18CBEnq48kppZxBLJCbC1VdL\nvf68eTKpR58+IrCffCJeczhJS5OnhU8+kXBW166Sp+jaVW44kydLH/6//hU2bJBe86NGyY2ooCC8\ntoWDOgm7MeYSY8wqY4zXGDMgVEYpilI7kpLgt7+V5GqfPpJAvOoqiSvHGj17ioBu3SrtAu68UwYZ\n/eUvkRHR9u3hd78TT/2tt6S1wogREipaswb+7/+khv622ySX0bGjxOLff18S2fWBunrsK4GfAp+G\nwBZFUepIcjLcc48IfM+eMkJzyhQRqlgjORluukk86Oeek/r9Ll0kIbxlS/jPb4wkTh9/XObCnTpV\nniY6doSbb5ZZpN54Q2wZOxb+8AfZdt99Mnl4LFMnYbfWrrPWbgBiOKesKI2PlBS4916JZwcmy77v\nPilLjDWMkRvQ9OnSOyYuTsJKl10myddI4HJJA7Q33pC5Ybt3h2uvlev20kvisS9aJCNdCwulp/zZ\nZ0tcvqQkMjbWBo2xK0oDJj0dHn5YEoY7dsiAnb/9TXqjxyKZmfDnP4uXPGyYNEc76yyp4/f5ImND\nq1YS1lq3Tq7V4sXyJHHVVXJjfOwxuZY33ggvvig233abhHZihZNWxRhj5gKVm2YawAL3Wmtn+/dZ\nAPzaWvvdCY5jp06devR9dnY22dnZwVuuKEqtWb5c4ssbN0oI4vzzo23RifF4JA7+6KMSC7/nHkkO\nu1yRtWP/fhHxp56SlsK33ir9cpo0kVzBCy9IOOmUU2RbXUYI5+TkkJOTc/T9tGnTolPuWFNh13JH\nRYkN5syB22+XqpAnnoCsrGhbdGKslYqWhx6Sro933y3livHxkbXD55Puk3/7m/SOv+46qa7p0EH6\n6Lzzjoj/+vVSznnDDXUfmBXtckeNsytKPWHcOAkdZGdLNcg994hHHKsYI/1gcnLEc37rLbkpPfmk\n1MhHCodDYvEffCA9aYqLJQE7caJMB3jJJWLjxx/D3r3Sn+byy2WAVCT92rqWO15kjNkOnAG8Z4z5\nMDRmKYoSbuLiZMTqihWQmysJwxkzIitAwXDmmTIhyNtvyyQdp5wiTx2RTmJmZcl5A2WbU6ZIUvX1\n1+VaPv205AoClUmnny7J1vLy8NumI08VRQHEA735ZgkdBDpK1geWL4f775cE8f33y0CkaHTA9Pmk\nhcEjj0i/nLvuEluaNpVt778vNfKbNkkYbMoUidefjGiHYhRFqccMHy5x45Ejxbt87LH6MSCnb18R\n1JkzxVvu0UNaB3i9kbXD4ZDWwQsXyvnnz5eb4/33SwfKCy6QMM3bb8N338m0iHfdJR5/yG0J/SEV\nRamvuN2SmFy0SDzMM86QySvqA2ecIQOMnn1WwiB9+0oyMxqBgmHDJA/w5Zcy21OPHlISuX271Oi/\n+qpcV4cDBgyAn/0stOWSKuyKovyIrl1FJG+9VRKtd98tc5LWB0aNEq/54Yelj87IkTKBeDTIypKw\n1po1UhrZr5/Uv2/ZIpU0jz4KmzfL52PHyiQloRiUpcKuKEqVGAO/+IV4khs2SHimvnjvxkgP9qVL\npSzywgtlgFE4wh41oVUrGXi1bh20bCle++TJ0vohJUUGRG3aJOJ+2WVSATR/fvBPGyrsiqKckIwM\nGWr/u9+J8Dz0UP2IvQM4ndIieP16eQoZMECePg4ejI49LVrI9fvhB+k7M3SotAleu1Y8+ltukW1X\nXy2J7H//O7jzaFWMoig1Zvt2GZRz8KD0UIm1yT1ORm6uJDM/+ACmTRPRdzqjZ8+hQzKg6fHHpT5+\n6lRJqoIkfz0eSEjQqhhFUcJI+/YyavXaa6WK5m9/i/2698q0ayfD/z/6SCpXhgyRRHG0SEmRwWEb\nNkgFzaBBMpI1N1duOMGOrFVhVxSlVhgjYYKvvpJ5TSdOlF4q9Yl+/WQ06B13iP3XXSclidEiNRUe\neEBCMsnJ0lP/178O3iYVdkVRgiIrS8r5OneWYfVffBFti2qHMRLLXrtWhLVnT3kCiWb+oEULSbKu\nWiVVSJ9/HtxxNMauKEqdee89iVffequEFqIZtw6WVavE/oIC+Oc/JUwTC+hk1oqiRI3cXCkpdDol\nft2mTbQtqj3WyoQfd90lzbseekhCI9FEWwooihI12rWT1rojRkgScOHCaFtUe4yRUaCrV0vlT69e\nMgK3vqEeu6IoIefDD2Vw0/33S222qadNvefNk5GigwfLZNsZGZG3QT12RVFignPPlcTqs89KaWQk\ne6aHknPOkZG3HTrI/Kf//nf9KO9Uj11RlLBRVCTtaVevlqZYXbpE26LgWbpU2gBkZsIzz0Quh6Ae\nu6IoMUVioiRSJ0+W4fOffBJti4Knf39p0NW/v9TBv/ZatC2qHvXYFUWJCDk5Umny4IMSt67PLF4s\nIaY+faQlQIsW4TuXeuyKosQs2dkyiOmxx+DOOyM/EUYoGTxYJsvIzJS+77NnR9uiY1GPXVGUiHLg\ngEz63KSJ1IxHu068rnz2mVQAjR4tlTNNm4b2+OqxK4oS86SnSxOutm2lkdi2bdG2qG6MGCHzrno8\n0mc9FnrWq7ArihJx3G4Ztj9pkiRVv/022hbVjeRkKYW8/34YM0Y892gGKDQUoyhKVHnnHbjhBqme\nGTs22tbUnY0bZfRqy5bS/bJly7odT0MxiqLUOy66SGrcr74aXnkl2tbUnVNOkSRx795SFjlvXuRt\nqJPHboz5M3ABUApsBCZZaw9Vs6967IqiVMvq1TKL0O23Sy/yhsAnn0hZ5KRJ0m89mK6XEe/uaIw5\nB5hvrfUZY/4XsNbae6rZV4VdUZQTsmMHjB8vIZlHHwVHA4gp7NkjoRmHA159tfahmYiHYqy186y1\nPv/bRUBmXY6nKErjJjNTJpdYskRmNaovk2afiFat4OOPpfZ94ECZeSrchPJ+OBn4MITHUxSlERIo\nh8zNhSuvhLKyaFtUd5xO+OMf4emnJacQ7qqZk4ZijDFzgVaVPwIscK+1drZ/n3uBAdbai09wHDt1\n6tSj77Ozs8nOzg7eckVRGjSlpXDFFbJ+800Z0NQQ2LxZBmh17QrPPw9JScduz8nJIScn5+j7adOm\nRX4GJWPML4ApwGhrbekJ9tMYu6IotaK8XEZ17twJs2bV/1GqAUpKZBq+xYul3PNEXS8jHmM3xowH\nfgNMOJGoK4qiBIPbDS+9BN26ycCfgoJoWxQaEhKkV/2NN8KwYaHvelnXqpgNQByQ7/9okbX25mr2\nVY9dUZSgsFYahy1cCHPnQlpatC0KHTk5EnK65x4p9Tx+timdzFpRlAZLQNy//FKqTBqSuG/ZAhde\nKL3e//EP8egD6MhTRVEaLMbA449Lb5lx42Sy6YZCp05ywyoslPbGu3fX7Xgq7Iqi1BuMgSeegCFD\nZBBTQxL3pk3h9ddlvtghQ6RjZLBoKEZRlHqHtRKPXrJEwjIpKdG2KLS89ppUzUyfDmPGaIxdUZRG\ngrVwyy3SY+bDD2V+1YbE11/LqNXOnVXYFUVpRPh80mRr3z54912Ii4u2RaFHq2IURWl0eDxw6aXg\ncknowuWKtkWhRatiFEVpdLhcMGOGDF664Qbx4hs7KuyKotR74uNlaP66dVLr3tiDAyrsiqI0CJo2\nhfffh88+k06KjZkGFo1SFKUxk5YmFTLDhkHbtjB5crQtig4q7IqiNChat5Z+7iNGSLngeedF26LI\no6EYRVEaHN26SfnjL34h9eCNDRV2RVEaJEOGwL/+Jc211q+PtjWRRYVdUZQGy/nnSyJ1/Pi6N9aq\nT6iwK4rSoLnuOrjmGplrtLg42tZEBh15qihKg8da+PnPZZTq9OngqEcurY48VRRFqQJjZOLoHTtg\n6tRoWxN+VNgVRWkUJCTA22/DK6/IPKoNGa1jVxSl0ZCRAe+9J7MUde4MZ50VbYvCg3rsiqI0Knr2\nhP/8RzpCbtoUbWvCgwq7oiiNjrFj4f77pcb9yJFoWxN6tCpGUZRGibVw/fUyb+rMmZJgjUW0KkZR\nFKWGGANPPw25ufA//xNta0JLnYTdGPOgMWa5MWapMeYjY0zrUBmmKIoSbuLj4c034e9/l6RqQ6FO\noRhjTJK19oj/9W1AT2vtL6vZV0MxiqLEJIsWwYQJ0su9e/doW3MsEQ/FBETdT1NAJ6VSFKXeccYZ\n8Kc/STL14MFoW1N36pw8NcY8BFwDFACjrLX51eynHruiKDHNzTdLs7A334ydZGowHvtJhd0YMxdo\nVfkjwAL3WmtnV9rvd0ATa+0D1RzHTq00ljc7O5vs7Oza2KooihJWSktl0NIVV8Bdd0XHhpycHHJy\nco6+nzZtWuiFvcYHMqY98IG1tnc129VjVxQl5tm6FQYPFq/9zDOjbU0UYuzGmK6V3l4ErKnL8RRF\nUaJNx44yQceVV8LevdG2JjjqWhXzBtANSZpuBW6y1u6qZl/12BVFqTfcd59Uy8yZA05n9OwIS4w9\nVKiwK4pSn/B6pfXAsGHwhz9Ezw4VdkVRlBCyZw8MHAjPPSfT60UDFXZFUZQQk5MDP/sZfPcdtI7C\n2HrtFaMoihJisrOlWdjVV4OvngzBVGFXFEU5Cb//PZSUwCOPRNuSmqGhGEVRlBqwbRucfjrMng1D\nhkTuvBqKURRFCRMdOsA//iH17bHeT0Y9dkVRlFrwy1/CgQMwfXpk+smox64oihJmHnsMVq+GF1+M\ntiXVox67oihKLVmxAs4+GxYvhs6dw3su9dgVRVEiQJ8+8NvfwrXXygjVWEOFXVEUJQjuukti7I89\nFm1LfoyGYhRFUYJkyxYYNAjmzYO+fcNzDg3FKIqiRJBOneDRR2VUaklJtK2pQD12RVGUOmAtXHIJ\ndOkSnpGp2gRMURQlCuzbJwnVGTNgxIjQHltDMYqiKFGgRQsZlTp5MhQWRtsa9dgVRVFCxs9/LiL/\nxBOhO6aGYhRFUaJIfj707g0zZ8Lw4aE5poZiFEVRokjz5vDUUzBpEhQXR88O9dgVRVFCzBVXQGam\nlELWFQ3FKIqixAB5eVIl89ZbMHRo3Y6loRhFUZQYoGVLePJJqZKJxsAlFXZFUZQwcOmlcNpp8Ic/\nRP7cIQnFGGN+DTwCtLDW7q9mHw3FKIrSqNi1S0IyCxZAr17BHSMqoRhjTCYwBtha12MpiqI0JNq0\ngYcegilTwOeL3HlDEYp5HPhNCI6jKIrS4JgyBZxOGZkaKeok7MaYCcB2a+3KENmjKIrSoHA44Jln\nYOpUyM2NzDldJ9vBGDMXaFX5I8AC9wH/jYRhKm+rlgceeODo6+zsbLKzs2tuqaIoSj2lZ0+4+Wa4\n/XZ4880T75uTk0NOTk6dzhd08tQY0wuYBxQhgp4J5AKDrbV7q9hfk6eKojRaSkqgXz94+GG48MKa\n/11UBygZYzYDA6y1B6rZrsKuKEqj5tNPZVKOVasgJaVmfxPtAUqWk4RiFEVRGjMjR8KYMRJvDyfa\nUkBRFCWC7NsnA5fmzpUa95MRbY9dURRFOQktWsCDD8Itt8i0euFAhV1RFCXCXH+9JFNffjk8x9dQ\njKIoShRYsgQmTIA1ayAtrfr9tG2voihKPeKmmyAuTjpBVocKu6IoSj0iP18GL330EfTvX/U+mjxV\nFEWpRzRvDn/8oyRSQ9kkTIVdURQlikyeLKL+4ouhO6aGYhRFUaLMkiXSZmDt2h+PSNVQjKIoSj1k\n0CAYOxb+539Cczz12BVFUWKAnTtlJOrXX8Mpp1R8rh67oihKPaVtW/j1r+E3IZi2SD12RVGUGKGk\nRMofn3sORo+Wz9RjVxRFqcckJMAjj8CvfgUeT/DHUWFXFEWJISZOlPr2Z58N/hgailEURYkxli+X\nKpm1a6FZMw3FKIqi1Hv69hXP/aOPgvt79dgVRVFiEJ8PHA5NniqKojQYHHVQZxV2RVGUBoYKu6Io\nSgNDhV1RFKWBocKuKIrSwKiTsBtjphpjdhhjvvMv40NlmKIoihIcofDYH7PWDvAvQVZdxg45OTnR\nNqFGqJ2hoz7YCGpnqKkvdgZDKIS9VvWVsU59+cdWO0NHfbAR1M5QU1/sDIZQCPutxphlxpjnjDGp\nITieoiiKUgdOKuzGmLnGmBWVlpX+9QXA00AXa20/YDfwWLgNVhRFUU5MyFoKGGM6ArOttX2q2a79\nBBRFUYKgti0FXHU5mTGmtbV2t//tRGBVqAxTFEVRgqNOwg782RjTD/ABW4Ab62yRoiiKUici1t1R\nURRFiQxhH3lqjBlvjFlrjFlvjPlduM8XLMaYLcaY5caYpcaYxdG2J4Ax5nljzB5jzIpKn6UbYz42\nxqwzxsyJhWqkauyMuQFsxphMY8x8Y8xqfyHA7f7PY+qaVmHnbf7PY+aaGmPijTFf+//PrDTGTPV/\nHmvXsjo7Y+ZaVsYY4/DbM8v/vtbXM6weuzHGAawHzgZ2AkuAK6y1a8N20iAxxmwCBlprD0TblsoY\nY84EjgAvBRLTxpiHgXxr7Z/9N8t0a+3dMWjnVOCwtTZmqqWMMa2B1tbaZcaYJOBb4EJgEjF0TU9g\n5+XE0DU1xiRaa4uMMU5gIXA7cDExdC1PYOe5xNC1DGCMuRMYCKRYaycE8/893B77YGCDtXartbYc\nmIH8OGMRQwz2zrHWfgEcf7O5EHjR//pF4KKIGlUF1dgJMTaAzVq721q7zP/6CLAGyCTGrmk1drbz\nb46Za2qtLfK/jEdydpYYu5ZQrZ0QQ9cS5EkN+AnwXKWPa309wy1k7YDtld7voOLHGWtYYK4xZokx\nZkq0jTkJGdbaPSACAGRE2Z4TEbMD2IwxnYB+wCKgVaxe00p2fu3/KGauqT9ssBQZxzLXWruEGLyW\n1dgJMXQt/TwO/IaKGw8EcT1jzkONIsOttQOQu+Ut/tBCfSFWM+AxO4DNH954A7jD7xEffw1j4ppW\nYWdMXVNrrc9a2x956hlsjDmNGLyWVdjZkxi7lsaY84A9/ie1Ez1JnPR6hlvYc4EOld5n+j+LOay1\nu/zrPOBtJIwUq+wxxrSCo7HYvVG2p0qstXmVJrp9FhgUTXsCGGNciFi+bK191/9xzF3TquyM1Wtq\nrT0E5ADjicFrGaCynTF4LYcDE/z5vunAaGPMy8Du2l7PcAv7EqCrMaajMSYOuAKYFeZz1hpjTKLf\nM8IY0xQYywkGW0UBw7F38FnAL/yvrwXePf4PosQxdvp/hAFOOIAtwrwAfG+t/Uulz2Lxmv7Izli6\npsaYFoHwhTGmCTAGyQXE1LWsxs61sXQtAay1/22t7WCt7YJo5Xxr7dXAbGp7Pa21YV2QO/g6YANw\nd7jPF6SNnYFlwFJgZSzZCbyKVBSVAtuQ6o10YJ7/un4MpMWonS8BK/zX9h0kVhhtO4cD3kr/3t/5\nf6PNYumansDOmLmmQG+/Xcv8Nt3r/zzWrmV1dsbMtazC5pHArGCvpw5QUhRFaWBo8lRRFKWBocKu\nKIrSwFBhVxRFaWCosCuKojQwVNgVRVEaGCrsiqIoDQwVdkVRlAaGCruiKEoD4/8B0t0oMUmZZ9gA\nAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x108aa0710>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure()\n",
    "plt.plot(x,t,'ro')\n",
    "for i in range(10):\n",
    "    prior_w_samp = np.random.multivariate_normal(prior_mean.flatten(),prior_covariance)\n",
    "    plt.plot(testx,np.dot(Xtest,prior_w_samp),'b')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compute the posterior. From the slides, the posterior covariance and mean are given by:\n",
    "\n",
    "$$ \\boldsymbol\\Sigma = \\left(\\frac{1}{\\sigma^2}\\mathbf{X}^T\\mathbf{X} + \\Sigma_0^{-1}\\right)^{-1}$$\n",
    "\n",
    "$$ \\boldsymbol\\mu = \\frac{1}{\\sigma^2}\\boldsymbol\\Sigma\\mathbf{X}^T\\mathbf{t} $$ (this expression assums the prior mean on $\\mathbf{w}$ is zero.\n",
    "\n",
    "Where, $\\Sigma$ is the posterior covariance for $\\mathbf{w}$, $\\boldsymbol\\mu$ is the mean, $\\sigma^2$ the noise variance, $\\boldsymbol\\Sigma_0$ the prior covariance."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "posterior_covariance = np.linalg.inv(np.linalg.inv(prior_covariance) + (1.0/sigma_sq)*np.dot(X.T,X))\n",
    "posterior_mean = (1.0/sigma_sq)*np.dot(np.dot(posterior_covariance,X.T),t)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sample some values from the posterior and plot them"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEACAYAAACj0I2EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXd4VGXWwH/vpCcEElpAOiiCFCmCIihRRFHXVbGuve2K\nBcuuvXyi6+piQ1Fsi7KiYkVU1gaWiKIiVXqXEnqANNJnzvfHyZAQKcnMJBnC+T3P+2Tmzr3vPVxm\nzj33vKc4EcEwDMM4+PHUtgCGYRhGaDCFbhiGUUcwhW4YhlFHMIVuGIZRRzCFbhiGUUcwhW4YhlFH\nOKBCd8695pzb4pybv5fP/uGc8znnGlaPeIZhGEZlqYyFPg44reJG51xLYDCwNtRCGYZhGFXngApd\nRH4Edu7lo1HAnSGXyDAMwwiIgHzozrk/A+tFZEGI5TEMwzACJLKqBzjn4oD7UHfL7s0hk8gwDMMI\niCordKAD0Bb4zTnngJbAbOdcXxHZWnFn55wVizEMwwgAEamSsVxZl4srHYjIQhFpJiLtRaQdkA70\n3JsyLydU2I+HHnqo1mUwOU1Gk9Pk9I9AqEzY4gTgJ6Cjc26dc+7qivoac7kYhmHUOgd0uYjIJQf4\nvH3oxDEMwzACxTJFS0lNTa1tESqFyRk6DgYZweQMNQeLnIHgAvXVVPoEzkl1n8MwDKOu4ZxDqmlR\n1DAMwwhzTKEbhmHUEUyhG4Zh1BFMoRuGYdQRTKEbhmHUEUyhG4Zh1BFMoRuGYYQRq1dDZmZgx5pC\nNwzDCCNuvRU+/TSwYy2xyDAMI0zYsgWOPBLS0yEx0RKLDMMwDlreeQfOPhvq1QvseFPohmEYYcIb\nb8AVVwR+vCl0wzCMMGD+fMjIgGBqh5lCNwzDCAPefBMuvxwiIgKfwxZFDcMwapmSEmjdGr79Fjp1\n0m1WbdEwDOMg5KuvVKH7lXmgmEI3DMOoZcaOheuuC34ec7kYhmHUIps3q2W+fj0kJsLEidC8OfTv\nXw0uF+fca865Lc65+eW2PeKc+805N9c596Vzrlkg/xDDMIxDnfHj4bzzVJl7vXDnneAJ0HdSmcPG\nAadV2PaEiBwtIj2Bz4CHAju9YRjGoYvInu6WTz6BlBQ47rjA5jugQheRH4GdFbbllnubAPgCO71h\nGMahyw8/QGRkmQIfNQr+/vfA54sM9EDn3KPAFUAmcFLgIhiGYRyavPaaWufOwa+/wrp1cO65gc9X\nqUVR51wbYLKIdN/LZ3cDcSIyYh/H2qKoYRhGBTIzoW1bWLECmjSBCy6A/v3httv080Di0AO20Msx\nAfgcGLGvHUaMKPsoNTWV1GByWw3DMOoAEybA4MGqzFesgClT0jjiiDTKqcsqU1kLvS1qoXcrfX+4\niKwsfT0cOEFELtzHsWahG4ZhlEMEuneHZ5+FQYPg+ut1MfSRR8r2qRYL3Tk3AUgFGjnn1qERLWc6\n544EvMBaYFhVTmoYhnEo8+OPUFwMJ5+scejvvw/LlgU/ryUWGYZh1DB/+Qv06we33AL33QdZWTBm\nzJ77BGKhm0I3DMOoQTZvhs6d4ffftbJiu3Ya4dK+/Z77WXEuwzCMMGfsWI1oSUqCV1+FU075ozIP\nFLPQDcMwaoiSErXIJ0/W3qEdOsDnn0OPHn/c1yx0wzCMMOZ//9MyuT16qKXep8/elXmghCIO3TAM\nw6gEY8bAjTdCQQGMHKm1W0KJWeiGYRg1wMKFsGiR+s/HjVPLvHfv0J7DLHTDMIwaYNQotc4BHn8c\nPvgg9OcwhW4YhlHNbN0KH32kKf7//S8cdRQce2zoz2MK3TAMo5p56SW48EKoX1+t87ffrp7zWNii\nYRhGNVJQoFUVv/tOx+TJ8MUXBz6utqotGoZhGPtgwgTo2RPatIFHH1WFXl2YhW4YhlFN+KsqPv00\nzJsHM2dWfjHULHTDMIww4quv9G+fPnDppTBtWvWezyx0wzCMamLgQPjrX7U0bnq6xp9XFrPQDcMw\nwoSfftIeoSefDLfeCrNnV/85zUI3DMOoBs46C844QzNEIyJg9OiqHW/10A3DMMKABQvg1FM1PHHw\nYFiyBBo3rtocptANwzDCgEsv1eiWH39UP/odd1R9DiufaxiGUcusXq3RLUcdpcW4hg+vuXMf0EJ3\nzr0G/AnYIiLdS7c9AZwFFAKrgKtFJHsfx4vX68XjsXuHYRh1n7/9Td0rX34J99yjKf+BUF0W+jjg\ntArbpgBdRKQHsAK4d38T3Na7N4vmzq2KXIZhGAcdv/8OEydCixYQG6ulcmuSSvnQnXNtgMl+C73C\nZ+cA54nI5fs4VrzAbT168Ozs2WapG4ZRZ7n2WrXO33pLFftxxwU+V2350K8B9ltqxgMMXL6cuWal\nG4ZRR1m5UjsQ7doFp50WnDIPlKASi5xz9wPFIjJhf/uNABYXFZH0yitccsklpKamBnNawzCMsOOf\n/4SLLoL33tPF0KqSlpZGWloa330HzZoFJkPALhfn3FXAX4GTRaRwP8eay8UwjDrNsmXQv79Gtlxw\nQeCRLRs3QrduMH8+tGxZfS4XVzr0jXNDgDuBP+9PmfsZ1nUA17/+uilzwzDqJI88AqecAllZcMMN\ngc/z6KNwzTW6qBoIB3S5OOcmAKlAI+fcOuAh4D4gGpjqnAP4RURu3Ncc2V2mcVSPKt1oDMMwDgrm\nzYOpUyE6Gt55ByIDdGSvXq3ummXLApelRjJFO3cW7rsPLrusWk9lGIZR45x2GhQWQvv28Prrgc9z\n5ZXQrh2MGKHvw7ba4oQJWs+gf38V2DAMoy7w9de6AOrzaRPoQFm4UOu+rFwZnDw14tTu0QPuvRcu\nuQSKi2vijIZhGNWLz6c1WjweGDUKGjYMfK477oAHHtAm0sFQY6uUt90GjRrB3XfX1BkNwzCqj3ff\nhYwMLcIVaHo/aN2XVatg2LDgZarRaos7dkDv3vDUU3DeedV6WsMwjGrD7zPftUtDDFu3Dmwer1cb\nSD/8MJx77p6fhW21xZ079W/DhtogddgwWLGiJs5sGIYRep55BvLy4F//ClyZA7zxBjRoAOecU7Yt\nL0+bSwdCjVjoKSnCwoVlBd5feknHL79AfHy1nt4wDCOkbNwIhx+ulvUPP6gPPRByc+HII2HSJOjb\nV7eJqHI/91y4+uowtdC3bIGuXWHzZn0/bJgulF59deB3IsMwjNrgb38D59SHHkyu5BNPwIknlilz\ngPHjYc0a+MtfApuzRiz0uDghPx+aNoVZs6BVKygogNRUOPNMePDBahXBMAwjJHz/vWaEjhmjij1Q\nVq1SRT5vnupD0IbSvXtrKOTRR4dxC7rGjYWsLA1ZTE7Wi9KtG2zaBMceqyE/tkhqGEY44/PBYYdB\ny5Ywc6Za6YHypz/BgAHaAMM/92mnwUknwX336bawTSzauVP/8c7p67594eOP9R8waRIMGaIrxj17\n1oQ0hmEYVeeOOzRS77ffglPmkydrUEj5RKSXXlKf+l13BSdjjfjQL74YEhPL/OUFBXDWWfD88/qI\n8eKL+n7t2pqQxjAMo2osXAjPPac6KyUl8HkKCuDWW3We6GjdtmSJpvu/8UbgdWD81IjLJT9feOop\neOEFDcnJydHPPB646ir4z39g9Gh45RWYPj24jCvDMIxQ4vOpm6VNG/j55+DmeuQRtfAnTtT3BQXq\ndh4+HK67bs99w9aHPmCAMGmSOvuHD1fXytSp/s81dGf6dHj8cfjpJ90vLq5axTIMw6gU112nLeW2\nbNGY8UDx10yfPVtvDqD6cMsWrbJY0Y0TtolFAwZoO6YePVSRr1qllRejotQNs3QpNG+uq8ft2mnI\njtV8MQyjtpk6FcaNU+9BMMrc59Mbw//9X5ky/+QT+N//4NVXg/PJl6dGFPrjj8P992vM5dq1Grq4\ncyf06VP2jysq0sVREX19xRWaFlsen8/H7NmzmT17Nj6fryZENwzjEGXjRhg6VPXWlVcGN9crr0BJ\nCdx0k75PT9ewx7ffhqSk4GX1U6O1XGbM0PZMl16qtQv+/W94+WUNBZo9u+yY+Hjo3FmTkV5/XX3t\ni+bO5ZVrriF1+XIA0jp25PrXX6eLhcYYhhFiiorgmGO06cSqVcEthKanq5s5LQ26dNE6MAMHakao\nP2xxb4StD738ObZtK8uCeustXeG94gpdCJ0/X+9WmZn6eb168Oc/w/jxPm4/pjfPzpu3+5HCh/Up\nNQyjehg2TPXTa69p4+dAEVEddswx8NBDuu2mm2DDBg3Z3p+rJWx96OVp0kTLRR5/vPrUs7J01feI\nI9SPnpmpn0VGalzmhAnQqlURRy4p2kNYDzBw+XLmzp1b0/8EwzDqMK+/rouUQ4YEp8xBbwpr1mg/\nCIA331S//BtvhM5vvgciUq1DT7F3pk8Xad9e5NprRbKzRcaOFUlMFPF4RLp3Fxk6VETvcT6BErmI\nCbKUjv6N8mF8vMyaNWuf8xuGYVSFb78VSU4WadpUZPv24OZas0akcWORuXP1/W+/6fv58yt3fKnu\nrJK+PaCF7px7zTm3xTk3v9y2851zC51zXudcr0BuJD6fj5iY2YwbNxevV+jZU630337TTNJFi2Dq\nVOGFF5bQuHER4OE9LqIXszmXiaRxAmlHdKSn+dANwwgBy5apRR4drZEtweTD+Hy6kHrHHeqJ2LZN\nfebPPadlT6qLyrhcxgGnVdi2ADgX+D6Qky6aO5fbevdm7YknknH6ABLn9eKWG1dzySVaX/h//4N7\n7lhPbo6X4Td3ZGjmHfRvdjsg5JHAZM7kzMiJfFMwnbfe8lBYGIgUhmEYSkaG1ldp104V7xlnBDff\nqFEapXfHHboIOnSorh1ecklo5N0nlTHjgTbA/L1s/w7odYBj93iM8Hq9MrxHD/GWuk0ExAsyvEcP\n2bnTK8OHi6Sk+GRw6/vkd1pIa9YI+OQYZsg1nQdIhw65pS4Yn9SvL3LMMSLNmok8/LDIxo1VeyQy\nDMPYtUukf3+RIUNEunYVycsLbr7589W1snq1iM8nctVVIueeK+L1Vm0eqsPlEmrmzp1L6vLle13g\nXLVqLqNHw8iRS1m5/kKu5G0+5Hzu4XFmcwzjl3zN3/62iX/+0xEV5cjO1pj2khL49Vc46igtDP/F\nF3+MYTcMw6hIcbG6WZKStILie+8Fl6Wem6vzPfWUWvtPPaUlct98M7ja6ZWlRqotjhgxYvfr5s2b\n0+QA+3ftmsfI2JPIyr+Qs/mEk/iOqQzibCZy990d6NNH49dvv139XTt2wGefaRPqZs00G2vYMLj2\nWrjmGq3DYBiGUR6fT3VEcTGsXw9PPqlGYaCIwI03am2WK6+E999Xn/nPP0NCwoGPT0tLIy0tLXAB\nVIjwcbl4S59Jyu+TQ4KM4P+kIRlydOO3ZOhQn3g8Gglz4YUirVqJnHmmrkpHROiUDRqI3HmnyLBh\numI9eLDIG29oJI1hGIbPJ3L77epqufxykUsv1W3BMHasSJcu6sL59luRJk1E5s0LfD4CcLlUKrHI\nOdcWmCwi3Sps/w64Q0Rm7+240n2k4jn8WZ8D/VmfRxzBsHHj9sj6rLjPZ237UXT0BD7/simpqfDN\nN5CbK0RH+2jatIi2bWM56yzH44/rY09REURFCYMGbadv3xzmzGnLDz84zjhDM1VPPVVryRwIn8+3\nO9a9Z8+eASUxhWIOwzBCgwg88IAGX1x2mbpDfvpJExkDZf58GDQIpk1T3TN4sLpvTjop8DkDSSyq\njHU+AdgIFALrgKuBc4D1QD6wCfhiP8fv9e7j9Xpl1qxZMmvWrN2WeWX22bxZre/6icXSIHqtRJAv\nkeRLhCuQegnFMnGiyLPPijRMLpJIt0vAK+CVxnGLZPTTK2TMGJF+/fTued11Ip9/LlJYuPc75MI5\nc2R4jx4yMT5eJsbHy/AePWThnDmVubmGdA7DMEKDzydy330i3bqJfPyxSEqKyKpVwc25c6dIx44i\nb74p8vvvIi1aiLz3XvCyEoCFHlTSUKVOsJ/EokDxer1yXdeBMoIHpRHbJIZ8iSJfoj3Z4pxPBgzw\nyeWdTpexXCXtWCVJ7BBHSWlkjE9uuklkzhyRp55S5Z6cLHLZZSKTJpWtcFfGNVQZOYOdwzCM0FBe\nmc+dq9FxU6YEN2dJiUbH3HKLyPr1mij5wguhkfeQUeizZs2SifHxIiCFRMl4LisNb/SKx+OT2NgS\n8VAsZ/GxzOAYmcTZciJpEkOORER4JSJCxDmRI48UGT1aM7peeEHkpJNE6tfXDNUHH/xdxsa1362I\nA8lOLS9noHMYhhE8Pp8+2XfrplZ0r14iTz4Z/Lx33CEyaJDIunUiRxyhRmKoCEShH/TO3GiKuZy3\nWEtb7oseSrNmhRQUePDh4TPO5BS+4UnuZBgv82jMKZx7bgaxsVC/vlZRu+027Wf62mvqW1+2TIvp\n/PJLfW7Nn0UvZvMA/+Qn+uE9+C+XYRxylJTAX/+qzemnTtWmEt27wz/+Edy848drga0xY7Tuy5VX\nBj9n0FT1DlDVQTW5XPbnypg+3StJ0askgmIBn0RSKK35XWIismTYMJ98+63Iyy+L9OwpEhcnEhmp\nkTERETp69hQZOdIrf+06QNIYIPfwmHRnnjRimxye9JWMGeOVZcvKVsX3tR5gLhfDqF3y8zWpZ/Bg\njXIbNkzk1FNFioqCm/fHHzV5KC1Nrf4HHgiNvOXhUHG5iJQtNn4YHy8fxsfLzUcfvcdi48I5c+Ts\nDtfLYSyWCAoFfOJxPklJUZ9527Yi//qXhhfde6/+58TGisTEqG8tNlYEfJIQtUkGR74qz8Z0l6s6\nD5Z/PbxGrrhCpGVLXfw468ztckqrB+Xl2I57XfQ8kJyGYVQP27eLpKaKXHCBSEGByL//LXL00SJZ\nWcHNu3ixhkmPH6+LoQ8/HHzI494IRKHXeD30UHKgcECfz8ecOXP59ddEXnvtCH77zeH1QqtWGrK4\nfr2WsGzSBM48UzuTfPstfPQRZGdDcjKkpPjYtMnLzp2RxMdD166Os8/WhITMTB83nPIYjdOPII2T\nSCSH45lOesv1jJp8L926eYiIsLBFw6hpVqzQ2ixnnqkJQ2+9pQmHP/0ELVoEPu/GjVre+6ab4KWX\n9G+o3SzLl0PTppCcfBA0uKhNfv1VL/5PP2mWWOfOcP75GkM6fboW6ImL02yxnj21CP3PP2uN9pQU\naNtWs8pWrICcHEhIKKFN/s9c4ZvMRbzLLhL5ieN5J+JEVrW4gKysWI49Vr8A/ftrBlliYm1fherF\nbl5GbZOWpun3jzwC11+v8eC3367GWqdOgc+blaVGX2qqGn333Qc33BAqqcvOccwx2rbzggtMoVeK\n9HS480744ANV7A0bwt13q3L//HOYOFHrL2Rl6f5Nm0JMjJYYyMlRpX/00RAZuZN109PZIa3JpT4x\nFNKO1RwWMYNzRvVn6NCOzJqlN5CffoI5c/SmcMwx0Lu3/u3RQ1vu1QWsTaBRm4joAuUjj2hjnFNO\n0UbM118PU6boQmig7NoFp5+uuuKnnzSl3995LVSIaFXGFi3ghRcOkhZ04URxMTz6KDzzjGaXRkRA\nr16q7M8+Wy3zTz/Vu/GMGZCfr9lku3ZBXh6A4PARgZc2rCWJTPKJY3VEc4qkISKOBg00iqZPHzj5\nZL05rFihRcVmz4bFi/Vzv5Lv3VtvFpWp/RBO+Hw+buttbQKN2iE3V5suL16sBlmHDtoZ7fLL1Ug7\n5pjA587LU/eN16stM4PNAN0Xjz2m+ub779WANIUeICIafvTAA+q/8vnUx96vn35Jhg6F2FhYu1ZT\ne/1jwwaol1BE5vY8irz1iKQILx7ExQAO57R4mMejoVNFRXq+yEj1z7dpo26f1q1VhnXrYOFCDZ08\n7DC1KLp1K/vboYPedGqL/blTZs+ezdoTT2So3ul2MzE+nrbTptG7d+8aldU4dFiyRJ+u+/aFF1/U\nJ+iPP9bf7scfq8szUAoKNIx5yxbYuVNvDl27hk52P598AjffrIZjSorqDI/HFHrQLFkC//633oVj\nY9XF4hx07KileW+8sWxRJSNDLe2ZM31MmZLNggX1yM+PQMSRkKD7JSRogfv0dLUikpN1voICtfiL\nivQGUpGoKL0Z+OvNFBToE0VSki7iHnaYVpFs107Pk5ioo149HeVf16sXfOnOA7lTTKEbNY2IKvAR\nI9TnfO21+tt66y1tLPHZZ/rEGyj5+drsYvFi/b1NmqR/Q82CBVoH5rPP1Md/8cW6BnDllabQQ8aO\nHfDOOzB2rEbDJCVpG6nsbFWQ3burT+2yy9Qv7mfTJpg7F778En78Ud0rubmqnFu21H2Tk1WR79gB\nmzersm/QQJOdfD79IuXk6KOeX9lHR+vNIT5elbyIKvnsbFX0iYlqmURFqfL2evVGsmuXjri4Pyr6\nyrxOSoKGDX28cMUZjFmcRjzaHqqiO8VcLkZNsnmzlsbetk0V+JFH6vaXX1Y36pQpwZXCzc7WAn7+\nRMNXXlEDL9Rs26bBEo8+CgMHqmunTx9dC4iONoVeLfz2m/YYfP99VbwtW6p1/vvvunAaHa1ukz59\ndCHmuOPg8MN1O6jinjhRH//8vni/pb5zp1aFbNQon/r1vSQn10PEkZ+vc2dk6PFxcUJERAk+n8Pn\ni6Cw0FFSonOI6JctLk7/+p8AcnKgeXP9YnfqpC6elBS18KOj9UbjHzk5kJPjY+XKreTleYiLa0Ju\nriMzE9LT88hYkUmuNCKKYhqTQWMy8Hp20Py03hx5ZCOaNgUpWcvMNx7jpA2/kuw28WvHZn+oomkY\nwSCi1RHvvFMXOx98UI0Yn09fv/uuZoO2bx/4ObZt09/wpk3q1771Vv1NhZr8fK3KeMIJapWfdZaG\nQd51l57PfOjVjM+nYYwTJ+qIjNS7amKi+tfnzy+Lbfd6dUW8Qwfo0kX9bu3aqVJ1TqNofvwRvp5a\nwKb0ItrLQpp41pGdXMLRpw3GSwobN6qffv16L8WFRSRJBjHsoiDGQ6NWrYiIjGPXLlXGu3b5SwaX\nyVqxa5P65fQz57QhSMuWKl+r5utZOvH/OHNjGg3YyvdHlrlU/O6Uc/PyyKUeGTRmG034JKYl0fc+\nTWJiOzZvplReYfXqQrZujSIy0kOLFo4WLdQt5H9CaddOf3CtW1euhHFtYOGX4ceqVdq4Zvt2+M9/\nytwp+flw9dX62/v4YzVYAmX5clXmPp/6y4Pxv++PkhL1+yckqFvnxhvVKr/wwrJ9TKHXICKqwKdM\n0fHLLxqCePzxmrgkoq6XuXNh9Wq1mP2RK343c+PGQvGORfTIT8dLJFnUZxPN2eRSaNc+ii5dHN27\n+5j39q2M/P0rCoknh0S20ZgxrXpyyrAH2L7dQ0aGWvLbtsHWrfo6P19Xyv2Lsj6ffomKi/Xvnn77\niv8/QgQlRHgKaN0+kQ4dhO1zX+XJrW8zgOlEIgd0p4joY+uGDTo2btQf3O+/61i9Wi2gZs1Uubdr\np6NjR32a6NgxuFZgwWDhl+FFXp62chs9WsOLb79djSnQxcpzzlFDYdy44NwiH3ygTZyPOkpj2ZOT\nQyH9HxHRG9OqVXpTeu89NRAr+vtNodcieXmanOSPOZ8xQ0MUjztOwxD9lnlGhv5HLl0KCxbks3EN\nxFJEA7KIJ49ISsgCaHk4Xm8s2dklFOwSfEQAguCIppgICohvFEdiYjSxsaq8y/8tb/n6C8n4rXav\nVxV7QYHKnZlZQN7GTHzEUUgMxUTh26MQmeOPSh+cE6KiPMTHq1XUvr0+ifTrpxEHzZuX/fD8lLd8\nu3btyYYNHlavLlPyy5frtVm1So/v1OmPo2nT6nkE9stnawHhgYgqu7vv1t/RyJF7rldNn66uimuu\ngYceCnzh3+eD666D//5X3Tgvvlh93y+Ahx/WUOgmTfTc772396cKU+hhhD9m9ZdfdBV7/nwdERGq\n9A4/HGJi0ol+5R5OLd5KDIXkUJ8NtOCbyDYUn/lXCgsbs3p1HjtWZJAvyeQTRywFePBRRBTFLpb4\neEf9+mpNNGyobpQGDVSxezxlPvbyw+fbU8FnZGwne+p3NPNCCVEUE0kJkWxwcbijehEZWY/CQsjP\nFzIziygoiMDrjcDrdVTtv1YAH5EUE00+nsh8Wh1en/Yd6pGSohEErVvrza9xY33K2LYNVq5UJb90\nqV5T53RRunt3vVl2765Wld+iD8ZdYtE6tY+IPvX+3//p0+Szz6qfufznzz6r0WjjxsEZZwR+rmXL\ndPFz82ZVrOecE7z8++Opp2DUKH19ySUanVPR6PETiEKvkSbRhyJ+xV0+ZlVE3QwLF6o1unJlCz6O\nv5pvsxqyig548NGcTeTE5jIgthHt2sHAgbF8/eLrDFs/n3rsIgIvXhxj2nXi6qee4/ffHStWqIXr\nj2PPyVG/fv36ZdEq8fGq8CIi9EfiH+qCaciamI548/IQ9PsjODLi4mkTn4CIWvwJCY7GjWN2/1vK\n/7tA5yssVKvfn3xVXOx37/gP8FBCDCXEQEkSS5bCkqUC7P976/HoF9//BLJypa5bfPyx/ylDrZxm\nKTkUrvuSE7O/o0vENP57ZBTDxpm75GBARFtLPvSQBguMGKF+5vL344wMjS9ft06fgstb7FWhpEQX\nPB99VI2I1auDq/FSGZ56Cv71L/0Ojx+vN5JQYxZ6LeP31564bDm5JDO1dV8G3T6S2MQj2LRJLYcV\ny7azYNpvxOVEkCcN2B7VCBfXnNzcSDweVXDl3S3+6BqfTxdK8/NV6fnDIOPiVMHHx6tfPyEB8GWz\nddkcmu3ajscJ2xLrcfixx5LUMLl0xb1sAH/Ytq/h33fz5q1kTP6Oht4YdlGPnSSRQRO20JBC6hET\n43bLXFysx/l/yBXj9P3zVvxaOSeIlJT6+D2lbirwuGKaNI2ieXO321ffsqWudfj/NmvG7kJq5nKp\nWbxeTax55hl9InvoIY3DrphE97//qTK/9FJVxDExgZ0vLU1j1tPT1V3z/PP7tpJDxYMPwhNPaE2n\nd99Vt+GBMJfLQUplXAR720ekzBdeWKh//cP/vrCwzL3i8+lC5ZYteqPYulVfb93qX1AVCguFBg1K\naNo0iqSGCyypAAAgAElEQVQkR1KSxqI3aKDD/9r/Y6rozik/yn++du06MkaNplOxo5AYCoilgFgW\nRdYjcsj5QGN27FBLOztbQzbz89WS2hvO6Q++vAtJnwIEh+DBBzi8lNcK+tuIiBCio73Ex/uoVy+K\n/HzHzp2q1Fu1gsYNt7Nl1nv02/ErzdwqVh0Bt44fbVZ+iMnK0sYyzz+v7rbbbtPkvYrKNTNTwxS/\n+Ub93CeeGNj51q2Dv/9dwxqdg7ff1mqM1UlJibpxvvhCM9FHjKi8f75aFLpz7jXgT8AWEeleui0Z\neA9oA6wBLhSRrH0cbwr9ICI3V7/4Gzeqe2jTpj1f+4eILlo2aaL+bv9o1Gjv75OSfPyj74Et34o3\nLq/Xw7p16jtfvFj96CtW6KLp1q3qCoqMVCuvqMiH1wsOH/5mXGWvnT+VGhH/dodfyTvniIlRWVu0\n0AikrKwstm2LZsOGOGJiHIcfrouy3brBUUf5gIU0aVJMr14W1lhZfD4tm/Hf/6pVfvrpqsj79v3j\nviKqdO+6S2srPfFEYNVKt29Xd8fLL6sh0qsXvP663sCrk9mz9YaRk6N1ZQYMqNrx1aXQBwC5wPhy\nCn0ksF1EnnDO3Q0ki8g9+zjeFHodJCdHFbs/ZNI/tm/f++vMTC03TMFmGhZvJ8FlklfPR/eB3Wnb\nvhENGkBBbjqzPxhLry3LSHBZLGmTxJVP3k+f47vQoMEfY9a9Xn1sXrlSx4oVwkf/+Z7I7MNYSxui\nKdLFZk80EpGo8fe+PGKlhCJiKKLUN7Xbf7/3305srFqQHTqoEsjcvp05368iZ1c7fERQP2ENJ57R\nioEnN6F3b12krYw74FCJdRfRtZ0PPtCkoMREjRu/5BJNdNsbCxfCLbfo9+allzSbsqpkZekC5Asv\n6EL777+rW+fKK6s3iiUrSyNZXn5Zvy8zZgQWG19tLhfnXBtgcjmFvhQYKCJbnHPNgDQR2WulYVPo\nBqjyzcyEnTt9zJixlNzcCBo3PoLsbA9ZWbBzp/DJixPokVFEFklkkkQWDVgV3ZSIhBZkZztiY/1u\nHyEychf16nlp3br+btdQgwaQl7OOX99+kV6bluMlgtkpR9LrgmF4acmMGZks+ymbTGlOPHlEU0Q+\nseQRS71ESEiIwustW9D1RwjtSdmGOPJpxibasIYNSQV0OuEMVq92rF4NXboIbdtuo3PnPM47rzVd\nu3r28AnX9Vh3r1cV2aRJOkpK1J1yxRWar7EvhbpunfrQP/9cXRQ33FB1//a2bRp6+MILGgW1bJmG\n0j7zjK6ZVBclJVoq5KGHdJ2qdWutzxJoD4SaVOg7RKRhuc/3eF/hWFPoxgE5ULhgr169yc2FX6cv\n5PXbRnDUmkzyqM+slE70veQG4hJbkZnJ7ptDeno2OTkRlJQkkJmpJQyKi4UEbwZNyCSWAiLw4iOC\nDOIoTG7Drl1R+HxqXftr5fjr2YjoDSknx7eHG+ePOKKifETJdlJK1hPtCtkY2RpvVDP694/ghBPg\n+ON9fPj3/oyZ/0udWXgVUbfYN9/o+P57XY845xxV5PtT4qBuvaefVlfMDTeoz7xBg6rJsHSpWuTv\nv6/Frnbs0GS2MWOqJ6LEj4j65f/xD/2+7NihCYZ+F09V2bhRb2IpKbUXtrhfjT1ixIjdr1NTU0lN\nTQ3RaY1DBecgIcHHJ/dezpvLyvnh18FtX37Bw3soQgf8URvk5wvD+1zEnYvWk136BLCDJF5p1Z3T\nb3mgtG4NrFmjP6rt23WBNju7TAaPx+HzeYikhFgKECCf2NKIGg/OCcXFPoppxO801l9GsY60NOHX\nXx1PPuljV04a01nKYKYwiO84kWkMXL6cuXPnVinWvTJtGEPt1hHR6zNnDsycqZ3AZs7UMNlBgzRC\nxe9uOBArVqhvfOJEtd4XLKhaRcOiIq0h/tprKs8ll8CQIXpTueceTamvjqJaoNfh2291oXPrVnUj\nPfecVnr8+9+r7tZ57bU0Ro1KY+lSTaQKUKhKNXpuA8wv934JkFL6uhmwZD/HimEcCK/XK8N79BBv\nuUAZL+g2r1dERGbNmiUT4+P/EFDzYXy8zJo1q1LnCaRpd0mJSHq6yOTJIv/8p1c6Jv1PejBbGrFN\nIigWh1fAK+ArHeXfeyWKfIkiRyIiSvypXRWG7h9NrnTosEtuvllkwgSRr7/2yrvvLpQvv5wnhYXe\nff5bJsbH77dB+b4+P9C/ecMGkRkzRCZO1AbLV14p0revSP36Ik2aiJx6qsj994t88onIxo2Vmnb3\n3J9+KnLGGSKNGok8+KDItm2VP15EZMECkdtvVzkGDhR58kmRa67RBvD33y+SmVm1+aqCzycydarI\ngAHaJPrNN0VGj9ZG85MnV20ur1evRb9+IgkJ2pz++uv1+0Z1NYl2zrVFXS7dSt+PBHaIyEhbFDVC\nhd+vPNDvVz7iiD2qNYYqizNYq7WinF+068PAv79KfklHJk/OYP5nm8iU1uSQiMO3Ox5e0TIKjmLi\nyceLo4B4IvHiJbKcK6csq9aDlyJiSEyE5s0jadwY6tUTVs/4mgFZ66lPDvXIpR45fNmiIRfeewcR\nEfDhY49w+fpVRFOCDw8FRPN2yw6cc+d9FBd7KCrSpw9d2ygbmzfraNhQfc4tWmhmc+fOOjp10uil\nqrJokZakHj9e5xw2DC64oHItGEXUev/wQx05OVq6ulcvte6/+UZdNbfcopFK1UF+vsr/3HMaLnz/\n/Rqlc/31GnX13ntwxBGVm2vrVi37++yzul5TUFDWcNofo15dUS4TgFSgEbAFeAj4GPgAaAWsRcMW\nM/dxvCl0o9LsT9mGU9LPvuQsL2MRMazkcBbTiZFNz6RBl6tYtMiRkQE+395/E5GRDo9HKCryEkUJ\ngocSysJ74uOhWTNHcnI2SfM+5DTvIqLwkkMiOSSyMDKRhHPOw+t1ZH76JS28UEwUEXiJpogNkT6S\nh55Fy5ZNiYrSBbvkZB1JSfq3WTN1e/gT1AK/Rlqc7vPP1a+dmanumMsuU5/6gdi1S5OApkzROYqL\nNXP0jDN0ofPVV1Wx33CDJhxVVwP2NWu0uuN//qMlsm+5RcvezpypLp7TT9ewyAO5dkpKNB597Fj4\n+mv1t4OW573hhj8WAwtEoVfJnA9kYC4XI4QE4jKpaQ4ko88nsn69yFdfeeXaazdI16450qCBT5yT\nvbhk1B3jKBIoEfCKxyMSHV0ikRQI+CSCYmnCZjmWn+TUiFfk8cdXynffzQnaPVVVfD6RlStFxo8X\nuewydYd06iRy660i06ape2F/ZGaKfPWVyIgRIiefLFKvnkhqqshjj4n8/LPIpEkiF12kLp9zztF9\nDzRnoOTkiIwbp+6cRo1Ehg8XWbpUP8vPF7nrLpGUFJEPP9z/PF6vyPTpZe6hFi1EGjQQOf54Pba4\neN/HUl0ul2AwC90INQdD/HagMm7aBKNHr2XKyIWsl2PIQP0Hst/qlyXUJ5d48tkZEY/P04DiYi15\n0Fy20JEVHMNMBvE1H3fLZsy8n/eZyFVZOYuKyqqGLlyoIYozZmhUx7HHqgU7ZMjea61IuZpGCxbo\nmDNH66n07q0RIv37q2vnxx/VOv/6a03o+stf1EqvDrdKdrZa0JMmacexE0+Eq67S5CB/tMrPP2u5\ngK5dNXpmbyn8Pp9WXP3wQ30yAc1qzs/X+a66qnJ9SS313zDqAOXdNiVEMZvefMdAXqw3lJ2+PkRH\nO6KiICfHS0FBWbZrGf73PhwlRJSGWXqJRMMqHQ0bwmHNcijZMJmTs6bQxrOQOe1acuG/HqNlu6NK\nO1ipktu+XaNa/BnDa9ZoYbTWrVXpdu6smZ7HHqs+96KisqSyjRt1//Jj+XKt09Otmyq2bt2gZ0/1\n2c+apUrz++81EWjwYHVpDBmimcmhRESjbL75BiZP1pvHCSdoqOXZZ++prLdsgXvvVUX/7LN7NqIA\njX2fMkU//+ILvQHExKiv/PTTNZnptNOq1tDFFLph1BH2tUDcqXtPFi5UhffDDzBtmgBFtGtXQFxc\nIitXeti0qaxblccjeDz6+9Nyx379sK/fpODxqNKPjdVQxNhYwbl8YmJ8NGmSQEKCIy6urLKmf2Rn\nq2LLyysr+dCsmVrp/tGmjWZNZmeXlXNYvFh97UVFmgDUr5+myR97bGg7Wvl8ejOZMQO++04VOWio\n5emn66hff89jCgs1Semxx9SyfvBB3WfnTrXCf/hBnx6WLdN1h8JCjUM/6yw47zy9EVVm0XdvmEI3\njDpEZdwhIqqkpk3zK3h9tO/WTRchFy/W/YqKysop79oF2dk+Ekq2kcwuHEIBMWTTgF3E45wH59we\nZY8jdtewL6BZ60TatI2nfXu1mqOjNREmKkqVV2Sk1gTKzCyLoNm0SZN80tM1Pvvww7X14VFH6ejW\nTRukhColf9eusvr5ixbpAuasWfoU0LevulNOOUWjUvZ2zpISjcZ55BF9ArnhBk1a+/lnvc6rVmmk\njsej/6ZOndQCP/VUdRkFu6AMptANw0DdIX7l/v332gKwQQN1gTRurAosP7+YIzMnM8i3kPpkk04r\nVtOeee5wtkUfSUKCh8aNhZ3rFtO8oIh84smlHnnEkeviEReDz+evXqnKPDZWFXpiotZoadGirIBb\nu3Zw5JGqtCtawVWlsFBvEjt26L8tPb1srF2rinzrVlXWnTvrDeOYYzRCZX81VURUUb/6qmasRkWp\nrGvWaARQYqKGF27dqjcjv6//lFOqx6dvCt0wjD+wZYsq+G++0QXG9HQtIRzr0mlZpMr8VKZwAe+R\n1j2d52dPZ+tWD1OnLmHh3x6lRVET1tNq91juWpMd0ZyUFA8NG6qVWlKiTwb+2vu7dqnbx9/TVkSf\nEjwerb9fv77+jY8v630bFaVDK2dqmGJxsb4uKiqz9ouLy0ItW7TYs659q1Z642jX7o/11H0+dQmV\nvwH4i7vNn6+Lsn6Zk5L0dXa29rf1d8bq108XbgN1o1QFU+iGYRyQzExNkBnz/C4WL4lBgAZk4I3y\nURKRwqmnRnDVVZCSMofNg0/YayJXy2+ncdhhvVm/Xm8Y/uGvsb9li7pZtmxRBe/vc+uvYe/vllVY\nqNvi4sqatPibr9Srp1ZxfLx+lpCg7xMSynzr/i5Z+fll9f/z88uUf1aWLu7m5upNxv8koc1MVLai\nIp2raVP13ffurS6hbt30b3U3v9gXptANw6gSeXk+/v3vdbz7blM2boyjSRPH1q2qKEWE5KiFPJj3\nEpfxNklkB5TIVVKiitXvU684du7UqBgtflamfPPyVEH72xh6vWWNWvyLvqDKOSKirE1hRITeGOrX\nZ3e/3eRkdbc0aqTzL1qkTy1JSZqUdNllgWW/Viem0A3DCJilS7X2+JtvqmuhY0dI+y6XBfOj8Ekk\niWynUf0F/PmKzgwechhHH60uj+qsLR4Kdu7UqJavv9Y2dgkJGoFy/vnqSglX+U2hG4YRNNnZ8MYb\nWk88Ph5uvNFHfv5KPvywKbNmNaBlS0f9+hq1UlKivuWuXXUR0j/atq0dV4XXqzcmfxXIX3/VKCD/\n4uXpp+si6cGAKXTDMEKGz6fJMqNHa8jfX/+qCTJTp6qyj47W9x07+jtGlY1NmzTmvF07teIPO0z/\n+kfjxhp5k5hY1gy8Mni96p7x98XduFHPt3x52bmbN9eIlvIj0IbStYkpdMMwqoXly1WxT5gAQ4dq\nVcD0dG3w7E+Hv+kmzR4FXaBcvVpD/jZs+OPYvl396nl5qtQbNNCFUY+nzCceEaFPAOWTl3Jz1e+d\nkqLjsMN04bJjx7Kng6o2xghXTKEbhlGtZGRo5uSYMdqE4a67NBt0zBh105x1Ftx9t8Z/VwZ/aGBW\nlt4EvN49hz9ZqXzkSyizR8MZU+iGYdQIeXmqwJ9+WsP97rxTsy9fflkt+f79tfZJnz61LenBSyAK\nPfzK1BmGEfbEx2s6/LJl2m7t8cc1c7JZM9120kkaRXLqqVo7xagZTKEbhhEwERGquGfMgFdegXff\n1eYVkZEa633++Tr+/GeYN6+2pa37mEI3DCNonIPUVI2Aefddjffu1EmzNn/7rSxk8MILNazQqB5M\noRuGEVKOOw4++ww++URbyHXpohmf8+ZpD9ATToCbb9YFViO0mEI3DKNa6N1bu/989ZXGsXfrprHt\nM2dqeGLnztqLs7CwtiWtOwSl0J1ztzrnFpSOW0IllGEYdYfu3eG997SU78KF2riiTRtNWpo2TTM3\nP/pIi3YZwRFw2KJzrgvwDtAHKAG+AIaJyOoK+1nYomEYu1m4EB56SBOS7rsPOnTQRKV27TRRqX37\n2pYwPKjpsMXOwAwRKRQRLzANGBrEfIZhHAJ07QoTJ+rC6ZdfwvXXw/DhGrvety88+qi5YQIlGIW+\nEDjBOZfsnIsHzgBahUYswzDqOr16qVJ//311ubz2mlrsM2dqFcTvvqttCQ8+Aq6HJiJLnXMjgalA\nLjAX8O5t3xEjRux+nZqaSmpqaqCnNQyjjnHccRru+P332oR52zY491wt/PWnP8HIkVrvpa6TlpZG\nWlpaUHOELPXfOfcvYL2IvFxhu/nQDcOoFCK6WHrPPZq0lJKiCUpjx2os+6FEjddycc41EZFtzrnW\nwJfAcSKSXWEfU+iGYVQJn08TlB54QDsJpadr4a+nngq+yfTBQm0o9GlAQ6AYuF1E0vayjyl0wzAC\noqhISwr885+qyAsKNASyf//alqz6sWqLhmHUSXJytLLj00+r9X7TTfDYY7XXwLkmMIVuGEadZutW\n9a+/+aZ2Jvr8cw2DrItY+VzDMOo0TZvC66/DkiXaxu7oo7U1XklJbUsWHphCNwzjoOPww2HOHBg/\nHt56SxX9pElWPsBcLoZhHNTs3AlDhsCCBVqL/fnntTDYwY65XAzDOORITta6MPffD4sXw+DBcNll\nsHZtbUtW85hCNwzjoMfjUYU+cSJER8PmzdCzpy6gZmXVtnQ1hyl0wzDqDIMGae31zEw4+WTYtAk6\ndlQ3TFFRbUtX/ZhCNwyjTtGypdZZdw5WrChride1a91fOLVFUcMw6iQ+Hzz8MLzxhrbD27wZ7rxT\nM06fflobbYQztihqGIZRisejCv3f/9bCXoWFMHcuXHMNDB0KF10Eq1cfeJ6DCVPohmHUaS6+WJtW\nDxsGr76qCn35cnXB9OkDf/877NhR21KGBlPohmHUefr2hR9+gGee0Zrr8fH6d9Ei2LULjjxSKzkW\nFNS2pMFhCt0wjEOCDh1g+nRte3fddVouoFkzreY4bZoq/E6dNPPU56ttaQPDFkUNwzikyM2F88/X\nSo3vvQcJCWWfTZumC6fFxfDkkxoGWVtYtUXDMIxKUFysRb2WLYMvvoCkpLLPROCDD+DeezWG/Ykn\noFu3mpfRolwMwzAqQVQUjBunvvWTT4aMjLLPnIMLL9SKjqefrhEy11yjXZPCHVPohmEckjgHzz4L\np50Gqakap16e6Gi45Ra14lNStFTv/feHdykBU+iGYRyyOKedjy6+GE48Edav/+M+SUnw+OMwbx5s\n3BjepQRMoRuGcUjjnDajvv56GDhw38lGrVqpm2bKFI1rP+oo9bWH0xJhsE2ibweuBXzAAuBqESmq\nsI8tihqGcVDw4oswciSkpUG7dvvfd+pUuOsuiInRYwYODK0sNboo6pw7DBgO9BKR7kAkcHGg8xmG\nYdQ2N96oSnrQIFi3bv/7Dh4Ms2fDzTfD1VerL37WrJqRc18E63KJABKcc5FAPLAxeJEMwzBqj5tu\nguHDValvPIBG83i0mcbSpXDOOXD22XDeedpoozYIWKGLyEbgaWAdsAHIFJGvQyWYYRhGbXH77XDt\ntarUt2w58P7R0XDDDVqut18/jZq58kr4/fdqF3UPIgM90DmXBJwNtAGygA+dc5eIyISK+44YMWL3\n69TUVFJTUwM9rWEYRo1wzz1aofGUU+C776Bx4wMfEx8Pd9yhSUujRsExx2gEzQMPQPPm+z82LS2N\ntLS0oGQOeFHUOXc+cJqI/LX0/eXAsSJyc4X9bFHUMIyDEhG47z746itdKK1fv2rHb9umC6bjxmly\n0h13aEx7ZajpTNF1wHHOuVjnnAMGAUuCmM8wDCOs8Mep9+unPvKqVmNs0kSrOM6fr8d27qzlejdt\nqh55g/Gh/wp8CMwFfgMc8GqI5DIMwwgLnIPRo6FRI7j8cvB6qz5HixaajLRwoVZy7NIFbr0VNmwI\nsaxWnMswDOPAFBZqbZdOnWDMGFX0gbJ5s1rur78Of/mL+utbtdpzHyvOZRiGUU3ExMDHH8Mvv8Aj\njwQ3V7NmqtCXLoV69aBHD63Rvnx5cPOaQjcMw6gk9etrud0334SXXw5+vqZNddF0+XK10AcMgAsu\n0E5KgWAuF8MwjCqyerUq3//8B848M3Tz7toFY8dCr15w4onW4MIwDKNG+PlnzQydOlVL64Ya86Eb\nhmHUEP36aeTKWWdVXxhiVQk4U9QwDONQ56KLYOVKVerff79nf9LawFwuhmEYQSACV10F2dkwcaIW\n7AoF5nIxDMOoYZyDV1+FHTvg7rtrVxZT6IZhGEESEwMffaTjrbdqTw5zuRiGYYSIBQvg5JO1mFev\nXsHNZS4XwzCMWqRbN21jN3SoVlqsacxCNwzDCDH33gszZmhD6cgAYwkDsdBNoRuGYYQYrxf+9Cct\n5DVqVGBzmMvFMAwjDIiIgAkTYPJkrftSU5iFbhiGUU0sXAgnnQTffqv+9apgFrphGEYY0bUrPP20\nVlDMyan+85mFbhiGUc1ce622oHvrrco3xjAL3TAMIwx5/nntKzp2bPWeJ2AL3TnXEXgPELSfaHvg\nQREZXWE/s9ANwzjkWboUTjgBvv66cuV2ay1s0TnnAdKBY0VkfYXPTKEbhmEAb78NDz8Ms2Zp96P9\nUZsul1OAVRWVuWEYhlHGpZdCaipcf71WaQw1oVLoFwHvhGguwzCMOstzz2k44/jxoZ87aJeLcy4K\n2AgcJSJ/qF5gLhfDMIw9mT8fBg2CX36BDh32vk8gLpdQdCw6HZi9N2XuZ8SIEbtfp6amkpqaGoLT\nGoZhHJx07w733w+XXQY//KD1XtLS0khLSwtq3lBY6O8AX4rIG/v43Cx0wzCMCvh8MGQIHH88lLN5\nd1PjUS7OuXhgLdBeRPaaB2UK3TAMY+9s3Kh10z/6SBV7eWo8ykVE8kSkyb6UuWEYhrFvDjsMXn5Z\nXS/Z2cHPZ6n/hmEYtczf/gZFRfDf/5Zts9R/wzCMg5BRo2D6dJg0Kbh5zEI3DMMIA376Cc47T0Ma\nmzSxjkWGYRgHNXfdBatWwYcfgsdjLhfDMIyDlkce0SJe778f2PFmoRuGYYQRy5dD06aQnGwuF8Mw\njDqBRbkYhmEcwphCNwzDqCOYQjcMw6gjmEI3DMOoI5hCNwzDqCOYQjcMw6gjmEI3DMOoI5hCNwzD\nqCOYQjcMw6gjmEI3DMOoI5hCNwzDqCOYQjcMw6gjBKXQnXMNnHMfOOeWOOcWOeeODZVghmEYRtUI\n1kJ/DvhcRDoDRwNLghepdkhLS6ttESqFyRk6DgYZweQMNQeLnIEQsEJ3ztUHThCRcQAiUiIiIehb\nXTscLP/JJmfoOBhkBJMz1BwscgZCMBZ6OyDDOTfOOTfHOfeqcy4uVIIZhmEYVSMYhR4J9ALGiEgv\nIA+4JyRSGYZhGFUm4I5FzrkU4GcRaV/6fgBwt4icVWE/a1dkGIYRAFXtWBQZxIm2OOfWO+c6ishy\nYBCwOFiBDMMwjMAIqqeoc+5oYCwQBawGrhaRrBDJZhiGYVSBam8SbRiGYdQM1ZYp6pwb4pxb6pxb\n7py7u7rOEyzOuTXOud+cc3Odc7/Wtjx+nHOvOee2OOfml9uW7Jyb4pxb5pz7yjnXoDZlLJVpb3I+\n5JxLL41+muOcG1KbMpbK1NI5921pAtwC59wtpdvD6pruRc7hpdvD5po652KcczNKfzMLnHMPlW4P\nt2u5LznD5lqWxznnKZXn09L3Vb6e1WKhO+c8gN+vvhGYCVwsIktDfrIgcc6tBnqLyM7alqU8pYvM\nucB4Eeleum0ksF1Enii9SSaLSK1GFu1DzoeAHBF5pjZlK49zrhnQTETmOefqAbOBs4GrCaNruh85\nLyKMrqlzLl5E8pxzEcB04BbgPMLoWu5HztMJo2vpxzl3O9AbqC8ifw7k915dFnpfYIWIrBWRYuBd\n9EsZjjjCsKaNiPwIVLzJnA28Ufr6DeCcGhVqL+xDTtDrGjaIyGYRmVf6OhfNam5JmF3TfcjZovTj\nsLmmIpJX+jIGDa4Qwuxawj7lhDC6lqBPZsAZ6Jqknypfz+pSZC2A9eXep1P2pQw3BJjqnJvpnPtr\nbQtzAJqKyBbQHz7QtJbl2R83O+fmOefG1vajd0Wcc22BHsAvQEq4XtNycs4o3RQ217TUPTAX2AxM\nFZGZhOG13IecEEbXspRRwJ2U3XAggOsZdpZpLdC/NDHqDOCmUhfCwUK4rmi/CLQXkR7oDylsHm1L\n3RgfAreWWsAVr2FYXNO9yBlW11REfCLSE33K6euc60IYXsu9yHkUYXYtnXNnAltKn8z29+RwwOtZ\nXQp9A9C63PuWpdvCDhHZVPp3GzAJdReFK1tKE7r8vtattSzPXhGRbVK2OPMfoE9tyuPHOReJKsk3\nReST0s1hd033Jme4XtPS+k1pwBDC8Fr6KS9nGF7L/sCfS9fz3gFOds69CWyu6vWsLoU+EzjcOdfG\nORcNXAx8Wk3nChjnXHypJYRzLgE4FVhYu1LtgWPPO/anwFWlr68EPql4QC2xh5ylXz4/Qwmfa/o6\nsFhEniu3LRyv6R/kDKdr6pxr7HdTOK3fNBj19YfVtdyHnEvD6VoCiMh9ItK6NOv+YuBbEbkcmExV\nr6eIVMtA79jLgBXAPdV1niBlbAfMA+YCC8JJTmACGiFUCKxDozGSga9Lr+sUIClM5RwPzC+9th+j\nvpe6J+sAAAB6SURBVMDalrM/4C33/z2n9DvaMJyu6X7kDJtrCnQrlWteqUz3l24Pt2u5LznD5lru\nReaBwKeBXk9LLDIMw6gj2KKoYRhGHcEUumEYRh3BFLphGEYdwRS6YRhGHcEUumEYRh3BFLphGEYd\nwRS6YRhGHcEUumEYRh3h/wHXvFrSL/W/UwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x108bffe10>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure()\n",
    "plt.plot(x,t,'ro')\n",
    "for i in range(10):\n",
    "    prior_w_samp = np.random.multivariate_normal(posterior_mean.flatten(),posterior_covariance)\n",
    "    plt.plot(testx,np.dot(Xtest,prior_w_samp),'b')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compute the predicted mean for a test point $\\mathbf{x}_*$ and its variance using the following expressions:\n",
    "\n",
    "$$ \\mu_t = \\mathbf{w}^T\\mathbf{x}_* $$\n",
    "$$ \\sigma^2_t = \\sigma^2 + \\mathbf{x}_*^T\\boldsymbol\\Sigma\\mathbf{x}_* $$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "pred_mean = np.dot(Xtest,posterior_mean)\n",
    "pred_var = np.zeros((len(testx),1))\n",
    "for i in range(len(testx)):\n",
    "    pv = sigma_sq + np.dot(np.dot(Xtest[i,:],posterior_covariance),Xtest[i,:].T)\n",
    "    pred_var[i] = pv"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plot the predicted mean and variance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAEACAYAAACnJV25AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAF7pJREFUeJzt3XuQXHWZxvHnnUlMMmGFIJio0SAKiq5kMCVlFRZpdbms\nyt1S3AURKBZKQVl1VwrLGqzVKnFrMW7VoiICwuKyCisXa0FQM2xhieQyw82EoFsJoiRcIshlCaT7\n2T/6dNIz09PTl9Pdp898P1Wn0tO388tvZp7znvecPhO2BQDIj4FeDwAAkC6CHQByhmAHgJwh2AEg\nZwh2AMgZgh0AcqbhYI+IpRHxi4h4MCLuj4jzkvtHIuLRiFifLEd3brgAgJlEo+exR8QSSUtsj0fE\nHpLWSTpO0kclPWv7ks4NEwDQqDmNPtH2Vklbk9vPRcQGSa9LHo4OjA0A0IKWeuwRsZ+kYUm/Tu46\nNyLGI+LyiNgzpbEBAFrQdLAnbZjrJX3G9nOSLpW0v+1hlSt6WjIA0EMN99glKSLmSPqJpFttf7PG\n48sk3WL74BqPcVEaAGiB7aba3c1W7FdI+k11qCcHVStOlPRAncFlfhkZGen5GBgnY2ScjLOytKLh\ng6cRcZikv5V0f0SMSbKkCyX9TUQMSypJ2izp7JZGAgBIRTNnxfxS0mCNh25LbzgAgHbxydNJCoVC\nr4fQEMaZnn4Yo8Q409Yv42xFUwdP21pRhLu1LgDIi4iQO3zwFACQcQQ7AOQMwQ4AOUOwA0DOEOwA\nkDMEOwDkDMEOADnT8CdPAQCdNzpaXiq3W8EHlAAgoyIkiQ8oAcCsRysGADIijTaMRCsGADIpdjVf\naMUAwKxHsANAztBjB4AMSKu/LtFjB4DMiQkddXrsADDr0YoBgB5Lsw0j0YoBgEyJKU0XWjEAMOvR\nigGAHqtuxaSBVgwAZAitGADAFLRiAKCH0j4jRqIVAwCZMbUNI9GKAQA03oqJiKWSrpa0WFJJ0ndt\n/2tELJL0n5KWSdos6SO2n+nAWAEgVzrRhpGaaMVExBJJS2yPR8QektZJOk7S6ZKesv31iPiCpEW2\nL6jxeloxADCN2m0YqaOtGNtbbY8nt5+TtEHSUpXD/fvJ074v6fhmBgAASFdLZ8VExH6ShiXdLWmx\n7W1SOfwj4tWpjQ4AcqpTbRiphbNikjbMqKR/sn1TRGy3vXfV40/ZflWN19GKAYAapm/DSK20Ypqq\n2CNijqTrJV1j+6bk7m0Rsdj2tqQP//h0r7/ooot23S4UCioUCs2sHgBmgdFkaV1TFXtEXC3pSduf\nrbrvYknbbV/MwVMAmNnkNsydd9Z7dvMVezNnxRwm6X8k3S/JyXKhpHsk/VDS6yVtUfl0x6drvJ5g\nB4BJ6rdhpI4Ge7sIdgCYqhPBzrViAKCLOnk2TAUVOwD0yMzVusS1YgAAtGIAoJvS/mtJtdCKAYAe\n6VQrhoodALqkGwdOJSp2AOiJxqp1iYOnAABaMQDQDd1qw0i0YgCg6xpvw0gcPAWADOpmtS5RsQNA\nVzVXrUtU7ACQMd2u1iUqdgDomuardYnTHQEAtGIAoFN60YaRaMUAQFe01oaROHgKABnSjSs51kLF\nDgBdQMUOAH2uV/11iYodADqu9WpdomIHgAzoZbUuUbEDQEe1V61LVOwA0GO9rtYlKnYA6Jj2q3WJ\nih0AeigL1bpExQ4AHZFOtS5RsQNAj2SlWpeo2AEgdelV61JHK/aI+J6kD0naZvvg5L4RSWdJejx5\n2oW2b2tmAADQ77JUrUvNXY/9SklH1bj/EtvvTBZCHcCsduedvR5Bk62YiFgm6ZZJFftztv+lgdfS\nigGQe+m2YaReHTw9NyJOlbRW0udsP5PCewJAX8haG0Zqv2LfV9KTth0RX5H0GttnTvNaj4yM7Pq6\nUCioUCi0M3YAyJR0qvXRZKn4ctMVe1vB3uhjyeO0YgDkyuRqvTP99c7/MetIlvIXEUuqHjtR0gNN\nvh8A9KXuhHprGq7YI+IHkgqSXiVpm6QRSe+VNCypJGmzpLNtb5vm9VTsAHIp/QOmE969s62YdhDs\nAPKiu9V651sxADCrZbkFU0HFDgAt6mwLZtdauAgYAHRSdcWeVbRiAKBBWfwwUi0EOwA0oB966xX0\n2AGgSd3pre9aGz12AEhbv7RgKmjFAEAd/dSCqSDYAWAa/RjqEsEOAA3pl1CXCHYAqKnf+urVCHYA\nmKRfWzAVBDsAVOn3UJe6HOylUqmbqwOApuQh1KUuB/v5K1bowbGxbq4SABqSl1CXuhzsq8bH9Z0z\nzqByB5Bp/RzqUpeDfUDSyk2bNEbVDiBD+vkMmFo4eApgVstTC6aiq9eKKUm688ADdcIhh3RztQAw\nxeiodNVV0ubN5UWStmzp3XjS1NWK/VPvOFRnX3GFBgbaW22pVNK6deu0bt06+vUAmlap0vfbr/z1\nli35CXWpy8H+7fvv1l++8xANDkp77SUVCs33sx4cG9P5K1Zoy+GHa8vhh3OmDYC25KH1MllXr8ce\nYU1dXfmOBQtKWrx4QGecEfrSlyY+o1QqaWxsTKVSSVefdZa+ee+9u7ZIJUnnDw9r1bp1be8JAMi3\n/my/NH899q4G++rV1lVXSePj0kMPSS++WN1GiQm3I6Q5c6R99tmh+c/fpVNf/BftMbBaS3fs0Mcm\njflHCxZox2WX6aCDDtIhhxwyJeArGwZJNR9vRFbeA0Dz+jPQKzIe7NXrKpVKOn/FCh0/vqeu1ic0\nrmE9piV6XHtLmqty0E8emzWgohboRe2jJ3WyrtNHdZG+ES/rmHnzNDgwoNEDD9TZV1yhtycHaB8c\nG9N3zjhDhU2bJGnK443IynsAaF7/n/XSfLDLdleW8qp2W7t2rW8YGrKlCcv1Q0P+9rc3+rTT7AMO\neN6v1GOerxcs7bRUTJbSlCVU9By95L31hF+78B5/61tFb99e9HnDwy5WvX9RKt9XLLoRxWLj71Es\nFr127VqvXbt2wmPNvEdeTDcXQDetXm2fdpq9cqW9bNmUuOmTRXazedvsC1pdmgn2tWvX1nzOaq30\nsbrSe2nMQ3rKA9qRBH790JdKHtTLXqDnvUR/9Fvilz7//Ee8aZP9wgv1fzAaGadtP7B+vc8bHvYN\nQ0O+YWjI5w0P+4H165t6j34xU2jXmwugGyYHev+Gep8FeyNV7HTPOXf5ct9zzz2+5pprfH0SmKu1\n0kfpVi/WY36Fnq0K+Z2ePvCnTmKEPWeOvWSJ/e532+ec86i/Ou+9/q329zP6C5dqhPJM/5dGg72R\nKrcblXC9dcwU2nnbO2HPo7/kK9D7MNjt3SFx/dCQrx8a8rnLl0+p7Oo9Z6YQqX58tVb6NF3pgzXm\nOXraEbWDfepSew9g4sah5AG95L30lFfoHn9OX/cP9WF/cd5xvvzyDV6/vuiPH/QhP6Z9/bwWuFQj\n7BqpcrtRCddbRyOhneZGrF3troM9j/5RCfTly+099ywvvQ/kPgh2Sd+TtE3SfVX3LZJ0u6SHJP1U\n0p51Xl/zG9JulTrTxqGRjcd0PyDz59sDA5Vw3zlpqd4wNBb+jT4W2ulBveS5A3/2295W8qmn2l/9\natHvee0/+990tq/Uab5OH/GPdYyP2/8cj44WvWaNfe+99oYN9u9+Zz/6qP344/bTT5fbTTt3Tvtf\nnjDP7e55NPKcRgKzG6Fcbx2N7nmksYFir6A1tcI8X4HenWB/j6ThScF+saR/TG5/QdLX6ry+Y9/g\nmX4xOv3LVywWfcKbzvTHdaWXa8yL9ZiH9Kyll2bYAFTfX73RmNw6qvXaWs9pd5lu47N7XQN6yUN6\nzvtqqw/WuE/X5T5l7hf92c8+4ksvtS+7rOgj3nChr9OHfaOO9a06ynfovT7pzaf57ruLHhsr+pS3\nHuffaj//Qa/xE3qVn9Ye/uTB7/LOnY3vvdT7njQSyjOto1sbqLQ2ct34HciKfFfnXQ52l8N52aRg\n3yhpcXJ7iaSNdV7b4W93bzW7Z1BtphBZvdr+4Aef8GI97IV61nO1wwPa6ai599BusE/eM6l3jKLe\nhqrZPZdGHq+1jmIyxpc8f/5Ov+lN9gknPO4z537W39WZvlYf83/peN+qo/zleUf6qqt+4/Hxok89\n6Bg/otf6KS3yC5rvnWqupZTGxiON92jkOf3eUqoO8oUL7XnzysvChb0O3PwG+/ZJj2+v89qOfvOz\noNWqqJ0DyWkelGwmaKo3YFdctmnCL97gYK2AnrjxiKoNx8QN1HQbl+oNWPVz6p8N1dpSewMTyQZk\naKjkpUv/z0cOfsff0tm+Wqf4Bp3g23SkvzLvfb722ge9cWPRp7/9CG/XK/2yBmvOZRobj2aOM9X7\nucnKgXt7apDPnVv+mZo7t9cBO3uD/ak6r/XIyMiuZfXq1R37wehH7R5I7uY42vkF70Yf/+c/L/qt\ni27ywRrzIj3ludpRc+MRVcv0x05qbUBqB3/rG49icjpuucW1LMZ94onbfMEFm33uK87yj3SS/1tH\ne1SHe41WeNX8Q3zLLff5iSfsu+5a5+sXTD8XabWUOlX1Tw7xwcHymWkDA5XjW7NxWW1ppGrpfrBv\nmNSK2VDntW3/EORdVqqmTq+jnTOd7Ob6362uY6a5qN54LE4+RNf6nkerex2NbmCKDu30HO3wXtru\n5VrvIwa/689/fosvvbTov3r9l/zvOtk36AT/RB/wT/V+n/jm0/2rXxU9Pm4/+GDRHz/oA/69XuMn\ntbef1UI/p/n+5MHv8o4dRdf6EakV2NMFWUR56X2gZnXpfLDvJ+n+qq8vlvSF5HbPDp6i/7RzplMa\n7YU09n7S2HhMfo/j33Smj/3Qk16+3J43r2QlexsT9ywm7znU24BUjkPMtOGo91mPRo97sPRdsEv6\ngaQ/Stoh6RFJpyenO/4sOd3xdkl71Xl9U780mN0aPZuknWDu9NlSjY6x3fdo5pTfi+Yd5bcsuskH\nHvC8Fy60BwZ2H+eYvi01tXVV78B9RPngZu8DMS9L88Hes4uAAe3qh6tlduuqoDM9Z7rHKxfjWzU+\nPu2lsBt5TsXoqPS1r5Wv4Lp9u/TyyzP//yov52/mTKePru4IIBsqVx5dWbny6AEH6Jwrr5xw5dFG\nntOKyuV0x8el3/5WevFFqVjc/fjAgDQ42NgGIr8IdgAtSGOvIG2V0L/rLumRR8qBXyrNxrAn2AHk\nWHWF/9BDu8PdznMrh2AHMItMbuXs3Fm+f84c6fnnezmyNBHsAGaxVavKQf/009LWrdL8+dIzz/R6\nVO0i2AFgQiVf+RunUr+GPMEOABOsWiXdeGP59vh4P4Y7wQ4AU1Qq+M2bd1fwW7b0bjzNIdgBYFr9\nGfAEOwDMqL/aMwQ7ADSkf6p3gh0AmjI6Wl4qt++8s4eDqan5YM/eVZMAoIsKhfJSqdqXLevhYFIy\np9cDAIBeKxTK/2a7cm8cFTsAqHbl3q/VOz12AJgkW313euwA0LZK9V65vXJlDwfTAoIdAGqYHO79\nhGAHgAb0U9VOsAPANPq1JUOwA0Ad/RjuBDsAzKDfwp1gB4AG9FO4E+wA0KB+OVOGYAeAFmW1aifY\nAaAJ/dCSIdgBoElZD3eCHQBakOVwT+UiYBGxWdIzkkqSXrZ9aI3ncBEwALkWTV2qq+F3bfoiYGld\nj70kqWD7Tym9HwD0heorQa5c2esrQZal1YqJFN8LAPpGFk+BTCuMLemOiFgTEWel9J4A0Hey0GtP\nq8f+GtuPRcS+ku6QdK7tuyY9xyMjI7u+LhQKKmRl8wYAKWqv1z6aLBVfbrrHnvpfUIqIEUnP2r5k\n0v0cPAWQW537q0s9+AtKETEUEXsktxdKOlLSA+2+LwD0kyyd/th2xR4Rb5T0Y5X77HMkXWv7azWe\nR8UOYFZJ5/TH5it2/pg1AKQo/ZYMwQ4AmdJ+1d67DygBABK9/tASHyoCgJT1+kAqrRgA6ILWWzK0\nYgAgM6pbMt1ExQ4AXUDFDgA50KuDqFTsANAlrVXtVOwAkCm9qNqp2AGgi5qv2qnYASBzul21U7ED\nQJc1V7VTsQNAJnWzaqdiB4AeaLxqp2IHgMzqVtVOxQ4APdJY1U7FDgCZ1o3rx1CxA0CPULEDQA50\no89OxQ4APTRz1U7FDgCZ1+mqnYodAHqsftXefMXO3zwFgJyhFQMAPdDJdgytGADIgOnbMRw8BYC+\n0amqnYodADKidtXOwVMAmPVSacVExNGSVqm8ofie7YvTeF8AyLtOtGPabsVExICkTZLeL+mPktZI\nOtn2xknPoxUDADOY2o7pzcHTQyU9bHtLeVBxnaTjJG2s+yoAgKT0r/iYRsV+kqSjbP9d8vUpkg61\n/elJz6NiB4AZpFGxc/AUAHImjVbMHyS9oerrpcl9U1x00UW7bhcKBRUKhRRWDwD9rboVs3z5qO69\nd7St90ujFTMo6SGVD54+JukeSR+zvWHS82jFAEADJrZjenDw1HYxIs6VdLt2n+64YYaXAQA6hE+e\nAkAGVLdjRkerz2dvvmIn2AEgg3a3YzgrBgBmPa7uCAAZkdblBWjFAEBGldsxtGIAYNajFQMAGZJG\nO4ZWDABkWAStGACY9Qh2AMgZgh0AcoZgB4CcIdgBIGcIdgDIGYIdAHKGYAeAnCHYASBnCHYAyBmC\nHQByhmAHgJwh2AEgZwh2AMgZgh0AcoZgB4CcIdgBIGcIdgDIGYIdAHKGYAeAnCHYASBn2gr2iBiJ\niEcjYn2yHJ3WwAAArUmjYr/E9juT5bYU3q+nRkdHez2EhjDO9PTDGCXGmbZ+GWcr0gj2SOE9MqNf\nvtmMMz39MEaJcaatX8bZijSC/dyIGI+IyyNizxTeDwDQhhmDPSLuiIj7qpb7k3+PkXSppP1tD0va\nKumSTg8YAFBf2E7njSKWSbrF9sHTPJ7OigBglrHdVMt7Tjsri4gltrcmX54o6YG0BgYAaE1bwS7p\n6xExLKkkabOks9seEQCgLam1YgAA2dDxT55GxNERsTEiNkXEFzq9vlZFxOaIuDcixiLinl6PpyIi\nvhcR2yLivqr7FkXE7RHxUET8NAtnI00zzsx9gC0ilkbELyLiweREgE8n92dqTmuM87zk/szMaUTM\ni4hfJ78z90fESHJ/1uZyunFmZi6rRcRAMp6bk6+bns+OVuwRMSBpk6T3S/qjpDWSTra9sWMrbVFE\n/K+kFbb/1OuxVIuI90h6TtLVlQPTEXGxpKdsfz3ZWC6yfUEGxzki6VnbmTlbKiKWSFpiezwi9pC0\nTtJxkk5Xhua0zjg/qgzNaUQM2X4hIgYl/VLSpyWdpAzNZZ1x/rUyNJcVEfH3klZIeqXtY1v5fe90\nxX6opIdtb7H9sqTrVP7hzKJQBq+dY/suSZM3NsdJ+n5y+/uSju/qoGqYZpxSxj7AZnur7fHk9nOS\nNkhaqozN6TTjfF3ycGbm1PYLyc15Kh+zszI2l9K045QyNJdSeU9N0gckXV51d9Pz2ekge52k31d9\n/ah2/3BmjSXdERFrIuKsXg9mBq+2vU0qB4CkV/d4PPVk9gNsEbGfpGFJd0tanNU5rRrnr5O7MjOn\nSdtgTOXPsdxhe40yOJfTjFPK0FwmviHpH7R7wyO1MJ+Zq1B76DDb71R5a/mppLXQL7J6BDyzH2BL\n2hvXS/pMUhFPnsNMzGmNcWZqTm2XbB+i8l7PoRHxdmVwLmuM823K2FxGxAclbUv21OrtScw4n50O\n9j9IekPV10uT+zLH9mPJv09I+rHKbaSs2hYRi6VdvdjHezyemmw/4d0Hcb4r6V29HE9FRMxROSyv\nsX1Tcnfm5rTWOLM6p7b/LGlU0tHK4FxWVI8zg3N5mKRjk+N9/yHpfRFxjaStzc5np4N9jaQ3R8Sy\niHiFpJMl3dzhdTYtIoaSykgRsVDSkarzYaseCE3cgt8s6RPJ7dMk3TT5BT0yYZzJD2FF3Q+wddkV\nkn5j+5tV92VxTqeMM0tzGhH7VNoXEbFA0hEqHwvI1FxOM86NWZpLSbJ9oe032N5f5az8he1TJd2i\nZufTdkcXlbfgD0l6WNIFnV5fi2N8o6RxSWOS7s/SOCX9QOUzinZIekTlszcWSfpZMq+3S9oro+O8\nWtJ9ydzeqHKvsNfjPExSser7vT75Gd07S3NaZ5yZmVNJ70jGNZ6M6YvJ/Vmby+nGmZm5rDHmlZJu\nbnU++YASAOQMB08BIGcIdgDIGYIdAHKGYAeAnCHYASBnCHYAyBmCHQByhmAHgJz5f/8NnEqGCB4k\nAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x108c63350>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure()\n",
    "plt.plot(x,t,'ro')\n",
    "plt.errorbar(testx,pred_mean,pred_var)\n",
    "plt.savefig('example.png')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
